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A Two-Dimensional Approximation to the 
Unsteady Aerodynamics of Rotary Wings 


ROBERT G. LOEWY* 


Vertol Aircraft Corporation 


SUMMARY 


A two-dimensionalized model is postulated for the representa 
tion of the aerodynamics of an oscillating rotary wing airfoil 
operating at low inflow; forward speed effects are neglected 
The resulting unsteady integral downwash equation leads to an 
equation for pressure distribution identical in form to that of 
classical fixed-wing flutter, but with a modified lift deficiency 
function. As a consequence, lift and moment expressions used 
in the fixed-wing flutter theory may be used directly for the ro 
tary wing airfoil in harmonic motion, except that Theodorsen’'s 
function of reduced frequency, C(k), is modified to include the 
effects of the number of blades in the rotor, the ratio of oscilla 
tory frequency to rotational frequency, and the inflow ratio. It 
is shown that the effects due to the vorticity shed by previous 
blades and/or in previous revolutions diminish as the product of 
reduced frequency and inflow ratio increases. Curves of aero 
dynamic damping coefficient for flapping and pitching motions 
are presented, showing that these effects (1) reduce the flap 
damping significantly at integer values of the ratio of oscillatory 
frequency to rotational frequency, and (2) make single degree of 
freedom pitch instability possible at certain values of this param- 
eter depending upon the location of the center of rotation. Ap- 
plication of these new coefficients to classical flutter analvsis and 


other rotor problems is discussed 


SYMBOLS 


pitch axis location (measured from midchord, positive 
aft) 

blade semichord (dimensional quantity ) 

F + 1G = Theodorsen’s lift deficiency function 

F’ + iG’ = modified lift deficiency function for rotors 

vertical distance between successive rows of vorticity = 
2ru/bQQ 

reduced frequency = bw/V = (w/2)(b/1 

> 


frequency ratio = w 


Presented at the Aerodynamics—II Session, Twenty-Fourth 
Annual Meeting, IAS, New York, January 23-26, 1956 

* Chief Dynamics Engineer. The work leading to this paper 
was performed while the author was employed at Cornell Aero 
nautical Laboratory 

The writer is indebted to Dr. T. Goodman for his valuable as- 
Sistance in overcoming some of the mathematical hurdles en- 
countered, T. T. Chang for his careful review of the manuscript, 


and 1D. Feign for his help in expediting the numerical calculations 


rotor revolution index 
= pressure 
= number of blades in the rotor 
= blade section radius 
= time 
= inflow velocity (steady, normal to rotor dis« 
= vertical component of induced flow (‘‘downwash 
= steady resultant velocity Qr 
= wake weighting function 
= ordinate in stream direction along the airfoil (positive 
aft) 
= displacement of airfoil perpendicular to rotor dis« 
= total circulation around the airfoil 
= line strength of continuously distributed vorticity 
= ordinate in stream direction (positive aft 
= oscillatory frequency 
= rotational speed 
= potential function 
= mass density of air 
= blade pitch (torsion) angle 
phase angle by which motion of gth blade leads that of 
the reference blade 
n,A,7= dummy variables 
Superscripts 
= dimensional quantities 
= time derivatives 
= functions independent of time 
Subscripts 
ad = quantity associated with reference airfoil 
n = rotor revolution index 
) = multibladed-rotor blade index 
l’ = refers to upper airfoil surface 
L. = refers to lower airfoil surface 


INTRODUCTION 


bes HELICOPTER ROTOR is subject to unsteady flow 
phenomena in many phases of its operation. Un- 
like fixed-wing aircraft, the rotor experiences oscillatory 
aerodynamic effects even in steady, gustless forward 
flight. This is fundamental for rotary wings since 
variations in relative airspeed occur at a blade section 
as the blade traverses the azimuth. In addition, flap- 
ping and blade pitch angles are made to vary with azi- 








82 JOURNAL OF THE 


AGRORAUTICAL 











Schematic elements of unsteady rotor flow field 


Fre. | 


muth in order to tilt the thrust vector in the desired 
Rotor 


blades, furthermore, have elastic, inertial, and aerody 


direction and eliminate overturning moments. 


namic properties in proportions which suggest suscepti 
bility to instabilities analogous to fixed-wing flutter; 
cases of classical rotor flutter have, in fact, been en 
countered, although it has not generally been a major 
problem largely because of the practice of mass-balance 
ing blades to reduce control forces. 

Still another problem should be considered in discuss 
ing unsteady effects because of its widespread impor- 
Azi- 


muthal variations in the steady downwash distribution 


tance--namely, vibratory rotor blade bending. 
over the rotor disc in forward flight have been shown! 
to induce blade bending vibrations whose associated 
fatigue stresses are often high enough to design the 
blades structurally. Similar downwash distributions 
are possible in hovering with certain multirotor con- 
figurations. Since this effect involves resonances be- 
tween the harmonics of steady air load distribution and 
the natural bending frequencies of the rotating blades, 
the magnitude of unsteady aerodynamic damping is 
most important in determining these stress levels. 
Oscillatory airfoil coefficients have become a stand- 
ard part of the fixed-wing aircraft engineer's equipment 
because of their importance in such problems as flutter, 
gust loads, and, under certain conditions, flight stabil- 
ity. Despite the indicated influence of analogous un 
steady rotor aerodynamic effects on rotor instabilities, 
blade stresses, and helicopter fuselage vibration levels, 
unsteady aerodynamic coefficients are not generally 
available for rotating wings. Isaacs* * evaluated the 
lift and moment on a two-dimensional airfoil in a flow 
with sinusoidally varying air speed, for application to 


helicopter problems. The complex nature of the wake 
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behind and below the rotor, however, was ignored jp 
order to make the analysis tractable; using typical 
helicopter parameters in a numerical example showed 
that the results of the theory gave small differences fron 
what is obtained by quasi-static means. Reference 4 
presents a closed-form solution to the same problem 
Considerable research into the rotor flutter phenome 
non has been conducted, prompted in some instances 
by specific cases of rotor instability, and more often by 
the possibility of gaining insight into the complexities 
of rotor dynamics and aerodynamics afforded by such 
studies. Asa result of both experimental and theoret; 
cal programs (for example, see references 5 and 6), it 
appears that more specialized treatment of rotor aero 
dynamic coefticients must be made than has been done 
heretofore. Better representation of the troublesome 
induced effects could proceed in two directions, which 
for convenience can be considered independently; (1 
accounting for three-dimensional effects in the plane of 
the rotor, such as the vorticity which results from the 
variation of circulation along the blade span and trails 
with circular form behind the blades, and (2) taking 
into consideration the influence of vorticity which has 
been shed and blown below the rotor disc and which 
must be passed over by succeeding blades and/or in 
succeeding revolutions. Cases of so-called ‘‘wake 
excited flutter,’’ reported in 1945 (see reference 7) and 
more recently encountered in tests at Cornell Aeronau 
tical Laboratory and the NACA (see references 6 and 8, 
respectively), strongly hint that the effect of vorticity 
shed by previous blades and in previous revolutions is of 
importance, at least for low inflows, since ‘‘wake 
excited flutter’? could not exist without them. 
Accordingly, the approach taken in this paper has 
been to account for the vorticity which lies below the 
After 


a brief description of the unsteady rotor flow field and 


plane of the rotor in a two-dimensional analysis. 


discussion of the two-dimensionalizing assumptions, a 
new set of unsteady aerodynamic coefficients is derived 
for application to rotor dynamics problems. It should 
be stated that this attempt is intended only as the first 
step in the development of a more satisfactory oscillat- 
ing airfoil theory for rotors. Components of airspeed 
in the plane of the rotor, for example, as experienced 
by a helicopter in forward flight, have been neglected. 
It is also appreciated that the range of helicopter, pro 
peller, arid fan parameters is such that the three-dimen 
sional effects mentioned above are undeniably signifi 
cant for many cases. The results of the theory pre 
sented herein are expected to be useful, however, in 
studying unsteady effects on rotors operating at low 
inflows as well as vibratory phenomenon at larger in 
flows in which the oscillatory frequency is significantly 
higher than the rotational speed. Limitations result 
ing from basic assumptions are discussed further in the 


text. 
(1) Rotor UNSTEADY FLOW FIELD 


A rotor in hovering or vertical flight trails a tip vortex 
which is blown axially downward so that, if otherwise 
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UNSTEADT 


undisturbed, 1t would form a contracting helix as shown 
in Fig. 1A. 
sumed to be constant (an assumption which is generally 


If the inflow over the rotor disc, u, is as- 
satisfactory for performance calculations, but not satis- 
factory for computing blade stresses), then the fluid 
that comes off the trailing edge of the blades makes a 
radial elements (see 


with horizontal 


Now, if there is an oscillation in blade ef- 


helical surface 
Fig. 1B 
fective angle of attack, blade lift will alternate also, 
and as a result of these changes in lift, vortices will be 
shed continuously at the blade trailing edge. These 
vortices fall along the horizontal radial elements of the 
helical surface shown in Fig. 1B, so long as the oscilla 
tions in angle of attack are small. Fig. 1C illustrates 
this helical sheet of shed vorticity. It should be noted 
that vorticity is considered to be on the helical surface 
shown in Fig. 1B; the vertical displacements from that 
surface shown in Fig. 1C represent the strength of the 
vorticity at a particular azimuthal and radial position. 
The variation in this vertical displacement (hence vor- 
tex strength) around the azimuth corresponds to the 
history of the motion of a given blade element at a 
fixed radius; the variation of shed vortex strength in 
the radial direction at any fixed azimuth angle is a func- 
tion of the variation with blade span of (1) blade chord, 

amplitude of the oscillation in effective angle of 
attack, 
Hemholtz theorem) vorticity cannot begin or end in 


space, a variation of shed vorticity in the radial direc- 


and (3) relative air velocity. Since (from the 


tion implies that there are trailing vortices at constant 
radii similar to and inboard of the tip vortex. These 
trailing vortices have been included in Fig. 1D. 

The schematic drawings of Figs. 1A through 1D in- 
dicate pictorially the complexities of attempting to ob- 
tain a complete representation of unsteady rotor aero 
One means of simplifying the problem is to 
“pitch” Con- 
sider, for example, the case of very high uw in relation to 


dynamics. 
consider extreme values for the helix. 
2, and the opposite, where uv is very low compared to 2. 
When 27u/QQ (the vertical spacing between adjacent 
helical surfaces of shed vorticity) is very large, then one 
would expect that all shed vorticity beyond a small frac- 
tion of a revolution would be too far below the blade in 
question to have an effect. Under these conditions, it 
would be sufficient to account for only the attached 
vortex sheet within that fraction of a revolution, as in 
Fig. 2A. On the other hand, when 27u/QQ is very 
small, all the sheets of shed vorticity tend to pile up on 
each other, and the effect of that vorticity close to the 
blade in question (shed by the several previous blades 
and/or in the several previous revolutions) is of more 
exists beyond a still 
the blade. 


importance than that which 
smaller azimuth angle on either side of 
This situation is depicted in Fig. 2B. 

The first condition is one which would exist at high 
rotor thrust coefficients; in this respect it may be of 
consequence in the phenomena of rotor stall flutter 
with high dise loadings. The second condition cer- 
tainly occurs at blade pitch angles near zero and is of 


priine importance in “‘wake-flutter."’ Vorticity shed 
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Schematic representation of unstead | 


by previous blades and in previous revolutions, as 
shown in Fig. 2B, may also be of substantial influence 
in improving other rotor calculations involving unsteady 
aerodynamics where the spacing between rows of vortic 
ity is not too large. This case that of low inflow is 
considered in this paper, and one of the results which 
definition of 


the analysis is intended to obtain is a 


“low inflow” in this regard —1.e., to determine under 
what conditions the vorticity beneath the rotor disc is 


important. 


THE AERODYNAMIC MODE! 


2) FORMULATION O} 


In arriving at a model which is mathematically tract 
able, the assumption has been made that for the case 
of low inflows, only the vorticity contained within a 
small double azimuth angle straddling the blade is of 
real consequence. In addition, it is assumed that the 
chord, amplitude of the oscillation in effective angle ot 
attack, and relative air speed vary slowly enough with 
span or in a compensatory way so that what occurs 
aerodynamically at one blade radius station is essen 
tially duplicated on either side of it. Asa result of the 
first assumption, the azimuthal angularity of shed vor 
ticity with respect to the blade may be ignored, as well 
as the horizontal distance that a spanwise element of 
shed vorticity is moved out of the vertical plane of a 
blade section at fixed radius as a result of azimuthal 
rotation. The second assumption is tantamount to 
saying that the flow problem at a given blade radius is 
two-dimensional. Thus, the portion of the circular 
cylindrical surface which is determined by (1) a par 
ticular blade radius, (2) the azimuth angle on either side 
of the blade section (within which the shed vorticity is 
of importance), and (3) the vertical distance spanned 
by a given number of rows of vorticity, can be regarded 
as a plane—one in which the two-dimensional unsteady 
aerodynamic problem may be attacked. This is shown 
3a, with an arrow indicating that the rows of 
From 


in Fig. 
vorticity under the rotor disc extend to infinity. 
the original postulation of small inflow, it is clear that 
the vertical angularity of the wake with respect to the 
plane of the rotor blade, 6, is negligibly small. Further- 
more, consistent with the idea that only the vorticity 
near the blade section has an important effect, one may 
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AIRFOIL "a" 
Vv=Ar’ 
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au , 
‘ 2rtu/an x VY = anuson 
a 2rtu/ON ee 2r7u/ON = 
27Tu/QN 2Ttu/QN 
a. b. 


Fic. 3. Two-dimensionalized models of unsteady rotor flow 








allow the planar rows of vorticity to extend to infinity 
in the horizontal direction in order to achieve mathe- 
matical simplification. These final modifications are 
shown in Fig. 3b, which is a two-dimensional model of 
unsteady rotor aerodynamics for a single-bladed rotor 


operating at low inflows. 


(3) DERIVATION OF THE UNSTEADY AERODYNAMIC 
COEFFICIENTS 


The general approach used here is the same as that 
of Schwarz.’ Only the effects resulting from oscilla- 
tory motion of a thin airfoil are considered, and the air- 
foil itself is represented by continuous vorticity, Y,, dis- 
tributed along a straight line. Fig. 4 is a repetition of 
the aerodynamic model presented in Fig. 3b, but show- 
ing the case of a multibladed rotor and the notation 
adopted for the mathematical analysis. 

Since not all the vorticity has been shed by one blade, 
the vortex rows must be identified by two indices; 1 
indicates which revolution a given row is associated 
with, while g reveals which blade has shed the vorticity. 

If the airfoil is sinusoidal, the downwash 
over the reference airfoil, v,, and its vorticity (either 


motion 


bound or shed) may be written, respectively, as 


ral =— a tol 
— = Cal 


Va = Val ug 
where the barred quantities are complex and independ- 
ent of ¢. The vorticity along the reference airfoil, y,, 
and that in the wake are differentiated from each other 
by writing wake vorticity as yoo for the attached sheet, 
and y,, for the sheet shed by the gth blade in the uth 
preceding rotor revolution.* 

The induced velocity at a point x’ on the airfoil, re- 
sulting from an element of vorticity of strength, y,,,, at 
a general point in the wake is given by the Biot-Savart 
Theorem as 


— &’)dé 


Tals 


dvq(x", t) Sa tlet!  HIY2 2(fp’\2¢ 
Qa) (x £")* + (nO + @)*(h')*; 


where ( is the total number of blades and h’ is the verti- 
cal distance between successive rows of vorticity. The 
total downwash over the airfoil can be expressed con- 
veniently by writing the integrals involving the bound 

* Note that 0 for the blade index is associated with the refer- 


ence airfoil (see Fig. 4). 
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vorticity and the row of vorticity attached to the refer 
ence blade separately from those for the rows of vortices 
below the plane of the rotor, as follows: 


eb 7] 7) . of 
Lf yalé’, tidé yool(é’, tdi 
OAS, B) _ + + 
Di | ae I ff — £! 
awWiLe A é Jb x 
j—% . as , 7 oS 
} ae Yngl&’, t)(x’ — &')dé 
— , -/\9 ) , By 
j=l wot (v’ — &)? + (nO + @)*(h 
. i / ,S 
.- vy nol’, t)(x’ — E dé 
ne , -/ , ) 
le (x — &)? + wh’)? ) 


The first two terms on the right-hand side of this in- 
tegral equation are those that arise in classical fixed- 
wing flutter; the third and fourth terms contribute the 
influence of all the vorticity below the plane of the rotor 
disc. 

It should be stated that this is a first-order theory, 
so that such effects as that of the wake upon itself have 
been neglected. The oscillatory flow induced at the 
airfoil in the free-stream direction, for example, will af- 
fect the airfoil lift and moment as well as the vorticity 
shed at the trailing edge. Since the shed vorticity is 
being accounted for as perturbations in the main flow, 
however, this component of induced velocity will be 
small compared to the free-stream velocity (except 
when the spacing approaches zero) and must react 
with small, first-order airfoil displacements. Because 
of the physical unrealities associated with the standard 
mathematical representation of the problem (such as 
zero thickness of wake and airfoil, zero displacement of 
the wake from horizontal planes, and infinite velocities 
at the core of a vortex), the accuracy of the results as 
h’ — Ois questionable. The effect of the horizontal 
component of induced flow, therefore, has been neg- 
lected in this analysis. 

It may also be noted that the effects of viscosity in 
damping the shed vorticity could be accounted for by 
considering only a finite number of rows beneath the 
reference airfoil, or by multiplying the vortex strength, 
y’, by an experimental “‘decay’’ function related to 
‘time’ through the rotational speed and the number o/ 
rotational cycles corresponding to a particular vortex 
row. This has been omitted since vortex decay time is 
long compared to the rotational periods of most rotors 

The vorticity shed by the gth blade » revolutions 
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q= BLADE NUMBER 
n=00;9=Q-/ n= SHED VORTEX ROW NUMBER 
Fic. 4. Aerodynamic model for a multibladed rotor showing 


notation for mathematical analysis 
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oo can be related to that of the reference revolution as 


a 


ic lows: 


therciore Y na Y 


Since a vortex is shed at the trailing edge of an airfoil 
with each change in total airfoil circulation, I’’, the 
following may be written for the reference airfoil 

F dV, ’(t dt ae 
yoo (b, thdx - —iwl,’e di 

at 
In addition, a vortex shed from the trailing edge of the 
ith airfoil is related to the change in circulation of that 
airfoil by 
je! dV,'(t)dt 
dé — 
dt 


2ar’g 
Ai. ee 


t Q 


—iwl,’e'dt 


If the amplitude of the circulation of the gth airfoil is 
equal to that of the reference airfoil but leading it in 


time by the phase angle y,, 


What is shed at the trailing edge of the gth airfoil at 
time, ¢, will be at a point in the wake described by 


t’,m some time, Af, later; where 


Af (29/2)(g/Q) + [(E’ — b)/Or’| + (22/2 
Thus, 
badly 
2ar'q : “ 
5 > = , t }dé —twl, ee dt 
O 
€ 
" ag - >) ‘O 2 < 2 oO - 
Yng\S )¢ dé 


‘! 
wer shai lai i (dt/dé') (2 
and for the special case where g = 0; 1 0 

Foole’ — iel'y’e~ ele — O18 de de" 3 
Substituting Eq. (3) into Eq. (2), yields 

Vaelt ) Jt’ e'*%e tw[(2ern/Q)+ (2 QO (4 


Eq. (4) is the phase relation between the strength of a 
vortex element at a general point in the wake and that 
in the plane of and shed by the reference airfoil. Re 


calling that d&’/dt Qr’, Eq. (2) may be written 


= 0, (2 x(w/2)(q i2n Q)n , 
(1wl cr je" %e é é 


Then defining the quantities, 


mAw/0: kA (w/2)(b/r') 2 m/r: FS (Ty /d)e 
one obtains 

Vng(é’) Ghee OE lg Seale, Se" (5) 
Substituting Eq. (5) into the integral downwash equa- 


tion above Eq. (1) gives [since v(x’, ¢) d,(x’)e | 
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eh , 
Yal€ dé 
(x ~ : — << Me a 
2r |. vr’ —¢& 
. , af oO l 
¢ ” dé te ~ ~ 
(7 eae _ ap BF went oree y 
. r' —€& q 0 
. mi (x! — &)dé’ : . 
| : - a 2 oO 
. r—é€& + (nV +q)"h 


[ _siies Mh 4... 3 
+ n°QOrh’? 


Str 


Nondimensionalizing 
the primes to indicate nondimensionalized quantities 


by dividing by 6, and dropping 


one obtains 


~7 Le r-—€ . t=? 
oO l 
ikI » 3 er e - V3 > x 
F é mS ( — E\dé ~ 
— — kT D0 e~?™™ x 
. (vx — &)*> + (nQ+q)h 


Str 


The last two integrals are of the form, 


[ e "(x — FS dé 

. (x — §)? +A 

and may be evaluated by the substitution a é 
—AX, and use of tabulated definite integrals (for ex 
ample, #11 of Table 103 in reference 10 Thus, 


> RE. : 
( (x — &)dé , 
i We 
(x — &)? +A 


Str 


v 2 , ie) ~— 0 \ = 
a ee.‘ me Di 
n=0 
ae ee 2 
I 
or 
o—} 
7 x —) ¥ ) 2 ) = 2 i) 
— oe . 
y=1 0 
Nv we iRx pd P 2 


Since kQh is always positive, the summations over ” are 


convergent geometric series which may be written as 


O-—1 
cite w/e 
1 
Do ime ' 
l—< s 
> pede a 
aw ime ( 
l1—¢ 7 


It should be noted that the closed form of the series 
the restriction that e “7? < 1, 
where kh = 0 are thus mathematically inadmissible; 
however, all singularities disappear in the final result, 
0 are in 


has Cases 


and numerical values corresponding to kh 
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cluded in the results for completeness. Eq. (6) may 


now be rewritten 


s] - & - 
—] Valé)dé 
24 ix-—¢ 


_— e Me mn k , 
kT - + skte*V 
J1 «x—& 
where (7) 


W = W (kA, m, Q, ¥,) 


Q-1 
( khQ i2aemi (0—g)/O.y 
: + zZ "es viii iets 
q 1 


AkO i2am 
ern =—% 


Ugl X) 


Note that the total effect of the infinite number of vor- 
tex sheets below the plane of the reference airfoil is 
contained in the third term on the right-hand side of 
this equation. The function, WV’, may be thought of 
as a weighting function for the vorticity shed by pre- 
ceding blades and/or in previous revolutions. 

The form of the integral downwash equation given 
in Eq. (7) may be solved directly by applying Séhngen’s 
inversion formula,!! as in the derivation of fixed-wing 
flutter coefficients. Séhngen shows the solution of an 
equation of the form 


“lt f(&)dt 
g(x) (1/27) 
ca 


x—€& 
to be 
. 2 fi-x f' 1+é g(t) |. 
Ft(s) — — dé 
T l+uJd-l 1 — &(x — &) 


if (1) is finite. Satisfaction of f/(1)-finite is tantamount 
to employing the Kutta condition in this case, where 


Yq at the trailing edge must be finite. Hence, 


ii-xf ef’ ji +¢ a, (edt 
Ya(x) = — 
T l+xid-! l1l—&x-—eé 


tkT f wre 5 
« l 


Zn ‘gs ~% “7 
oe & l+te™ 
a | dé (S) 
} uh we oe oe 


Performing the second integration on the right-hand 
side of the above equation, and then evaluating the cir- 
culation over the airfoil, , from 


+1 
r = e" / Valx)dx 
1 


yieldsf 


t A more complete derivation is presented in Cornell Aeronau- 
tical Laboratory Report No. 75. The fact that a solution for the 
unsteady pressure distribution can be effected using published 
results for an airfoil entering a generalized gust was suggested by 
Dr. N. Rott from the floor during the January, 1956, Annual IAS 
Meeting. This latter method, following the original approach of 
Sears (in the Journal of the Aeronautical Sciences, Vol. 8, No. 3, 
p. 104, January, 1941), was developed completely by Dr. N 
Kemp in subsequent correspondence with the author and provides 


an alternate means of development 
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P = | 
2 | i+é or 
—1 % }-¢ °° 
j | 2 
imk 15 [FT,?(k) + tHo'?)(k) |] + 
Q) 
where J,(k) and H,°)(k) are Bessel functions of the 
first kind and Hankel functions, respectively 
Bernoulli's equation for unsteady motion using non- 
dimensional lengths is 
(O¢/Ot) + (V’/b)(0g/Ox) + (p/p) f(t) 
where ¢ is the velocity potential and p the pressure. 
Since the airfoil has been replaced by a vortex sheet, 
the difference between the velocity on the upper and 
lower surfaces is y,(x, ¢). Thus, for any point on the 
airfoil, 


(Ogy/Ox) — (Og, /Ox) bya(x, ft) 


Px 
ov — OL b } val&, tid 
J -1 


Combining the above with Bernoulli's equation gives 
bu — pr = Ap = 


oe i ] 
_ Vya(x, t) + 6 valé, tidé 
o| Ya or Yals QS 





For sinusoidal oscillatory motion, then 


Ap(x*) odes iad 
V —7F¥,(x*) — ik Yal&*)dé* (10 
p J -1 


This equation can be evaluated in terms of downwash 





on the airfoil, @,(&), by substituting Eqs. (8) and (9) in 


the appropriate form. It can be shown thus that the 


equation for the unsteady pressure distribution be- 


comes 


— Ap(x*) 
pr’Q 
(22/m) [10(R) + 2J0(k) W] 
IT (k) + tH (Rk) + 2[Si(k) + iJo(k) |W 


f h—x* lise... 
2 ' 1—x* |[1 + 
T 4 i+ Wi - x*—f 


| M(E\dE (11 
where 


7 ] l—x*§+ V1-&2 V1 - xv 
2:(x*, £) = In 
») sk F -9? *o 

2 l—x% — V1 —-# V1 — x* 


Eq. (11) may be recognized as being in the same form 


x 


vt 


Str 


tRO(x*, & 


as the corresponding equation for the two-dimensional 
fixed-wing oscillating in incompressible flow. If the 


factor multiplying the first integral is written as 


(2/mr) [1 — C’(k, m, h)] 


the | 
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the function C’(k, m, h) is directly analogous to Theo- 
dorsen’s well-known lift deficiency function, C(k), given 
in re‘erence 12. Since the integrals required to solve 
for lift and moment remain unchanged, all the previous 
derivations of two-dimensional lift and moment on os- 
cillating airfoils (for example, see reference 13) apply 
for the case considered here. 

The function C’(k, m, h) is obtained from the defini- 


tion, 


at Cc’) 


D4 [Fo (Rk) + 2Jo(k) W(RA, m) | 
nr {Hi?(Rk) + tH? (Rk) | + 2[Si(R) + tJo(R) |W (Rh, m) 


Tt 


so that 
C'(k, m, h) 
HI? (k) + 2J,\(k) W( kh, m) 
Hy?(k) + tH (Rk) + 2[Si(R) + tJo(k) |W (RA, m 


(12) 


This modified deficiency function can, of course, be 
separated into real and imaginary parts, as is done in 


the classical flutter case, 


C’(k, m, h) F'(k, m, h) + 1G'(k, m, h) 


and the oscillatory lift and moment for a helicopter 
rotor at low inflows can be calculated by the following 
familiar expressions, using the notation of reference 13, 
where F’ and G’, functions of k, m, and h (or of k, m, h, 
Q, and y, in the case of multibladed rotors), have been 
substituted for and G, Theodorsen’s functions of k. 


Ly l — 72( 1/k)(F’ + jG’) 
a (1/2) — j(1/k)}1 + 2(F’ + 9G’); - 
2(1/k)?*(F’ + jG’) - 
o>) 
My = 1/2 . 


It will be noted that the moment expressions are un- 
affected by the wake since in reference 13 the moments 


are taken about the quarter chord. 


(4) THE EQUIVALENT SINGLE-BLADED ROTOR 


A number of general observations are pertinent re- 
garding the new lift-deficiency function given by Eq. 
12). It is clear that C’(k, h, m, Q, ¥,) approaches the 
Theodorsen function, C(k), as IW approaches zero. 
From Eq. (7) it may be seen that this limit, which ex- 
presses the negligible effect of vorticity below the rotor 
plane, occurs as would be expected, as the spacing be- 
tween the rows, /, becomes infinite. It should be 
further noted from Eq. (7) that, although the param- 
eters k and m are related through the expression 
k = m/r, the actual value of m is important in the ex- 
pression for k, whereas when it appears as a phase angle, 
only the noninteger portion is of consequence. 

For the case of a single-bladed rotor, one may write 


W = 1/(e*e**™ — 1) (14) 
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Now, the following facts are pertinent regarding a par 
ticular multi-bladed rotor whose blades are moving with 
prescribed relative phase, y,. The parameter & ex 
presses the number of oscillations in vortex strength 
which occur within a semichord length behind the refer 
ence airfoil. Thus, the value of & is directly tied to the 
motion of the reference airfoil and (since it is assumed 
that all blades are oscillating at the same frequency) is 
uninfluenced by the number of blades in the rotor. Be 
cause of the form of the aerodynamic model, however, 
the reference airfoil cannot differentiate between vor 
tices shed by preceding blades and those which it shed 
in previous revolutions. 

Thus, the fact that there are two additional param 
wake —namely, 


eters associated with the unattached 


h and m_ ~and that Wis a complex number (with two in 
dependent components) suggests that for a multibladed 
rotor under any condition of inflow and with arbitrary 
frequency ratio, m, there is a_ single-bladed rotor 
with equivalent values of h and m that has the same 
value of IV. 


Consider a two-bladed rotor: 


W=2 (Lt eee ee) / (Cre l 
For the special case where y¥ 0 
Woer = (14+ OO™) [Cr 


which can be equated to the function II’, for a single 


bladed rotor under equivalent conditions, thus, 


Wo-2 Wo-=1 lig e. = 7 
y =0 
igre l 
so that h, = ho, which can be written as u, us/2 and 
mM, m2/2. Likewise, if W 7, then 
h. = he or u, U>/2 
and 
m1, (m2/2) + (1/2 
Likewise, for a three-bladed rotor, 
2kh 1(4/3 ” ty hi ty 
W ] — ¢ ¢ ¢ 7 ¢ t 
V=3 h 2am 
t or 28 ae 
and, for y; Yo 0, 
W o<s 
y 0 
2Qki { ” 2 » 
¢ ¢ - € ¢ “T ] 
[ hs ,t(2 rm” is 1 | [e7**e! 13 ma Pais et (2 7 ” 1 | 
] 
Ce ee 
so that 
h, he (or u. = U;/3 
and Mm, = M3/3 


In general, for ‘‘collective pitch’? type disturbances 


i.e., Where y, = 0, 
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versus k as a function of inflow parameter (for m = 0 


Fic.5, F'’ 


(15) 


i collective khQ i2an 
€ e 


which always results in 4, = Ag and m, mg/Q. For 


“cyclic pitch” type disturbances, y, 22(g/Q)(w/Q), 


O-!1 
hO i2em kh 
i+e"er" 2 
, q 1 . 
W eyctic khO i2am (16) 
" illieg) — | 
The latter expression for ¥, expresses the fact that each 


blade is experiencing the same disturbance as it passes 
a fixed point in space. 

Since the unsteady air forces on a blade of an arbi- 
trary rotor can be expressed in terms of those of a single 
bladed rotor through the equivalent IJ’ function, as 
discussed above, the remaining general discussion 
shall be concerned with the case of a single-bladed rotor. 

For this case, the real and imaginary parts of the lift 


deficiency function are given by 
JiaA mf Y; + BJ,)B 


F’(k, m, h) 
(k, m, h A? + B 


Cith mp = Ut BWA + AaB 
a (R, , ) 
m,h A + B 


Jia + Jy, — 3A 
y; —_ Jie — Iva 


where A 
—B 


h kh 


€ al 
and a bh 


e — 


2 cos 2am + e 


2 sin 27m 


kh 
e — 


2 cos 27m + e 


Here again it is clear that as the spacing between rows 
of vorticity, 4, becomes infinite, F’ —~ F and G’ — G. 
On the other hand, since & and # appear as a product in 
the same terms resulting from the additional vorticity, 
C’(k, m, hh) > C(k) ~ 0.5 ask— @. Thisis physically 
interpretable as due to the fact that as the number of 
oscillations in vorticity within a given length becomes 
infinite, the distribution of vorticity below the rotor 
plane approaches symmetry about any point on the air 
foil, so that the contributions of this vorticity to the net 
downwash at any point on the airfoil approaches zero. 
Other limiting cases of interest occur when & and/or 


h approach zero. When h is zero, all the vorticity lies 
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in the plane of the airfoil; if, in addition, m is integer, 
J\(k 


C’(k, integer, 0 
a a iJi(k 


Ji(R) 
If hh — 0, m ¥ integer, and k — O (by virtue of r — 
then C’(0, Finally, if k — O and 


m — () or integer, 


m, 0) — C(O) — |. 


C’(0, integer, h) ~h/(h + x 


When / = 0, the meaning of any of these results, physi 


cally speaking, is not clear. 
thickness for both airfoil and wakes, and that the wakes 


The assumptions of zero 


are undisplaced from their trailing planes are all, oj 
course, unrealistic. 

Fig. 5 shows the real part of the lift deficiency func 
tion, F’, plotted versus the reduced frequency for vari 
ous values of the spacing and zero phase between suc 
cessive rows of vorticity. Beyond a certain value oj 
k, the function looks very much like a damped sinusoid 
which oscillates about the curve for infinite spacing 
The 


curves depends on the inflow ratio, or spacing, and the 


i.e., Theodorsen’s function. “damping” of these 


curve for h 0 oscillates undiminished between 0 and 
+1. Fig. 6 is a similar plot of the imaginary part oi 
the lift deficiency function, and displays similar charac 
teristics. Beyond a certain value of k, —G’ oscillates 
about Theodorsen's function —G like a damped sinu 
soid, apparently with the same period as the oscillation 
in F’, but between the limits +0.5. 

Both F’ oscillate with the same ‘‘period’ 
beyond the critical value of & for all values of the phase 
however, the 


and —G’ 


parameter, m,; as m changes, critical 
value of k and the psuedo-sinusoids are shifted along 
the R-axis. 


values of k of practical interest are in the lower range, 


Because of this periodicity and since the 


Figs. 7 through 14 are plotted for values of & from 0 to 
2.0. 
markedly the lift deficiency function in general; it is of 


Both phase and inflow ratio are seen to affect 


interest that below about k 0.5 and for m near 0.75 


the real part of Theodorsen’s function represents / 
rather well for all values of 4. When the spacing be 
tween successive rows of vorticity exceeds about 5 


chord lengths (/ 10), C’ does not differ very much 

















from C except near k = 0, m = 0. 
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Fic. 15. Pitch-damping coefficient versus frequency ratio as a 
function of inflow parameter 


For use in dynamic rotary wing calculations, tables 
of F’ and G’ as well as the real and imaginary parts of 
the corresponding values of LZ, and L, are given in 
reference 14. These results are presented for all com- 
binations of the parameter values listed below: 


/ 0, 1/4, 1/2, i, 2, 6, 30 

k 0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.8, 1.0, 1.2. 
1.4, 1.6, 1.8, 2.0, 2.5. 3.0, 3.5, 4.0, 4.5, 5.0, 5.5, 
6.0, 6.5, 7.0, 7.5, 8:0, 8.5, 9.0, 9.5, 10.0 


Note that values for / © (fixed-wing flutter case 
and # = O are included for completeness. 


5) APPLICATIONS TO ROTOR ANALYSIS PROBLEMS 


While the analytical attack employed in this develop- 
ment is two-dimensional, the coefficients derived herein 
may be used as a “‘strip theory”’ as well as in ‘typical 
section”’ analyses. For example, if the rotor geometry 
is given, and inflow as a function of radius and the rota- 
tional speed are known as operating conditions, then 
the parameter / is known for as many radial stations 
as desired. If, further, the oscillatory frequency, , 
is assumed, the remaining parameters m and & are de- 
fined. It will be noted that m is constant with blade 
span, but & varies as the semichord divided by the 
radius. 
may vary markedly over the blade span, and also that 


It is thus apparent that the values of k and / 


k, h, and m will have widely different values depending 
on the operating conditions and the phenomenon being 
investigated. 

From review of the basic simplifying assumptions, 
it may be concluded that the new coefficients provide a 
good approximation (1) for low inflows, where devia- 
tions from fixed wing theory are greatest, and (2 
where k > 1. The latter stipulation is analogous to 
cases where finite span corrections are unnecessary in 
fixed-wing flutter theory. It is clear from the propor- 
tions of the usual helicopter and propeller blades, how- 
ever, that, for k > | near the blade tips, w/2 must be a 


large number. 


Any of the standard methods of flutter analysis 
be used to solve the equations of motion arrived at using 
the present coefficients. Their application appears t 
involve the least error for the phenomenon first ri 
ported in reference 7 as ““‘wake-excited flutter’? whic! 
occurs at low inflows with w to 2 ratios as high as 20. 
Che oscillatory frequency in these cases seemed nearly 
equal to that of the fundamental torsion mode; on 
might have therefore expected a forced excitation 
resonant frequency where there is low damping, or 
single degree of freedom instability. The former js 
implied in reference 7, whereas the tests reported 
reference 6 revealed a region of instability. Figs 
through 19 were plotted to investigate the pt ssibility 
of pitch instability. They show curves of the coei 
ficient for damping in pitch for particular blade stations 
and pitch axis locations and various values of inflow 
The form of this coefficient is the same as that given i 


reference 15—namely, 


1/k)} (1/2) — a] — (1/2) + a](2G’/k 


where a is the pitch axis location measured aft of the 
mid-chord in semichords. The damping moment per 
unit span is obtained from the coefficient for an airfoil 
oscillating in harmonic motion with frequency w, through 
an angle 6, by multiplying by rpb‘w°@. These curves re 
veal that negative damping can exist when the oscilla 
tory frequencies are certain close-to-integer multiples 
The possibility of single degree of 


he 


of the rotor speed. 
freedom flutter is, therefore, shown to exist using t 
new coefficients. 

Fig. 17 is presented to emphasize the fact that, with 
the rotational axis on the quarter chord, the wake can 
have no influence on the airfoil moments. It is thus 
clear from these figures that as the pitching axis is 
moved aft from the leading edge, the ‘‘peak”’ values of 
negative damping become more stable, until there is n 
possibility of instability with the pitch axis on_ the 
quarter chord. As the pitch axis moves aft of the 
quarter chord to the midchord, negative damping ap 


pears again. 
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6. Pitch-damping coefficient versus frequency ratio as 
function of inflow parameter 
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Since “‘wake-excited flutter’’ has been encountered 
with blades which have pitch axes essentially on the 
quarter chord, it appears that such cases cannot be in- 
terpreted merely as single degree of freedom flutter in 
pitch. A more likely explanation is suggested in Figs. 
20 and 21, which show that flap-damping can go to very 
low (but positive) values at integer values of the fre- 
quency ratio. It thus seems likely that ‘“‘wake excited 
flutter’ is actually a classical type of flutter in which at 
least some of the aerodynamic derivatives are greatly 
modified due to the vorticity shed by preceding blades 
and/or in previous revolutions. This large change in 
the unsteady aerodynamic coefficients with frequency 
ratio at low inflow suggests itself as the mechanism by 
which it is possible to obtain two distinct types of classi- 
cal flutter with one blade configuration which has only 
two degrees of freedom, as reported in reference 6. 

Damped amplification factors have been applied to 
the problem of predicting rotor blade stresses by Goland 
and Flax" and others, using quasi-static aerodynamics 
to determine the value of damping. Fig. 22 is a com- 
parison of the amplification factor curves computed for 
an existing helicopter blade using quasi-static aerody- 
namics and the coefficients of the present theory for a 
single-bladed rotor. The quasi-static curves are those 
calculated by Daughaday and Kline and presented in 
reference 1. This example is typical of good blade de- 
sign, since the first three natural bending frequencies 
are removed from integer orders of rotor speed. The 
significant fact illustrated by these curves is that the 
aerodynamic damping, and hence the peak amplifica- 
tion factor, is critically dependent upon how close the 
frequency ratio, m = w/Q, is to an integer number. 
From Figs. 20 and 21, it is clear that the peaks of the 
curves obtained with the present theory could be higher 
still if m were closer to an integer number. The inflow 
ratio chosen for this calculation is roughly that of 
presently operational helicopters in hovering flight at 
maximum gross weight and was assumed constant over 
the blade radius. The effect of inflow can be to raise 
or lower these peaks, depending on the frequency ratio. 
In using the coefficients for a single bladed rotor di- 
rectly, it is tacitly assumed that only one blade has been 
disturbed. Actually, the blade bending moment prob- 
lem, being closely associated with the azimuthal dis- 
tribution of downwash over the rotor disc, falls into the 
classification of ‘‘cyclic’’ disturbances described in 
Eq. (16). It is intuitively obvious and can be easily 
shown, using Eq. (16), that the spacing, 4, between rows 
of vorticity for a Q-bladed rotor oscillating at integer 
multiples of rotor speed is effectively 1/Q that of a 
one-bladed rotor with the same value of inflow. From 
this it is clear that, when m is near integer, the more 
blades in a rotor producing a given force, the higher 
will be the amplification factor peaks. It should be 
noted that this effect is in addition to the reduction in 
flap damping resulting from the decrease in blade chord 
which often accompanies an increase in the number of 


blades in a rotor and which is compensated for, since the 
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4 


generalized exciting force is also proportional to the 


chord. 


(6) CONCLUSIONS 


Expressions for two-dimensional oscillatory airfoil 
coefficients have been derived for rotary wings, taking 
into account the shed vorticity beneath the plane of the 
rotor disc. The influence of these modifications to the 
fixed-wing theory appears only in a modified lift de 
ficiency function. The effects of inflow ratio, number 
of blades in the rotor, phasing between the motion of 
adjacent blades, and ratio of oscillatory to rotational 
frequency have been included and appear as parameters 
along with the basic ‘‘reduced frequency”’ of fixed-wing 
unsteady aerodynamics. The results all approach 
those for fixed wings as the product of inflow ratio and 
reduced frequency increases. For values of the fre- 
quency ratio whose noninteger parts are between, 
roughly, 1/10 and 9/10, the effects of the out of plane 
vorticity are significant until the spacing between ad 


When 


the frequency ratio is closer to integer values, the pres- 


jacent rows exceeds about five chord lengths. 


ent theory indicates the differences between rotary- 
wing coefficients and those for fixed wings to be sub 
stantial to somewhat higher values of inflow. 

Attention has been drawn to the fact that the lift and 
moment on a blade in a rotor with an arbitrary number 
of blades (whose motions are in arbitrary phase) can be 
obtained from the coefficients for an equivalent single- 
bladed rotor. 

From limited calculations, it 1s concluded that single 
degree of freedom flutter is possible in pitching motions 
for rotors operating a low inflow, when the pitch axis 
location and frequency ratio take on certain values. 
Similar calculations for flap-damping coefficients re- 
veal that the aerodynamic damping of rotor blades in 
motion perpendicular to the plane of rotation falls off 
abruptly as the frequency ratio approaches integer 
values. Furthermore, in general, the lower the value 
of inflow, the lower the minimum value of both pitch- 
damping and flap-damping coefficients. 

From these considerations, it is concluded that the 
quasi-static aerodynamic theory is inadequate for rotor 
dynamics problems such as calculating critical flutter 
speeds and/or damped amplification factors. The 
current trend toward increasing the number of blades 
per rotor, magnifies the importance of these effects, 
since an increase in the number of blades acts in general 
as a reduction in inflow. 

Correlative calculations using these coefficients are 
needed to determine their accuracy and the limits of 
their practical usefulness; data are available from rotor 
flutter and response tests which can be used for this 
purpose. Furthermore, this first and simplified ap- 
proach to the problem should be extended to include 
some of the important three-dimensional effects and an 


attempt made to account for forward speed. 


(Continued on page 144) 
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3 flow configurations \ vacuum spectrograph-monochromator 
atio and mera with a 1,470 A. radiation source was adapted for instal INTRODUCTION 
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of plane | formed by flow over a sphere, a cone-cylinder, wedges, and a graph, schlieren interferometer) are not feasible 
veen ad- vlinder with axis aligned parallel to the direction of the flow for low gas densities (pressures < 1 mm. Hg There 
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ition as indicated by the photographs is presented. Limita absorption of the ultraviolet is a function of the density 
» elt eel tions of the method are discussed of the oxygen. Therefore, the intensity of radiation 
number passing through a region of high density will be less 
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ore ( = constant in linear emulsion response equation : 
: ene ya ration. 
D = density of photographic negative _ 
it single E = calibration factor (In 10)/Tp,Cl Oxygen has a radiation absorption peak in the ultra 
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. . " ‘ ; aa — an If intensity Deam O! proper wave ieng which 18 
ents re- | = intensity of incident radiation ” tex ” I ERE ag : 
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other wave lengths in it. Since ultraviolet light is ab 





falls off I; intensity of radiation transmitted through developed 
integer negative sorbed by most materials, even in thin films, the optical 
a K = constant in microphotometer-recorder response equa- 
ie value . 
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h pitch- K, = constant in microphotometer-recorder response equa 
tion 
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t = pressure 
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R = outside radius of the nonuniform axially symmetric 
| flutter region of flow, see Fig. Al 
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PHOTOGRAPHIC PLATE 


Fic. 2. Caleium fluoride optical system schematic 


elements must be made of material which absorbs very 
little of the radiation, yet has the other qualities neces 
sary for making good lenses. A vacuum system is 
necessary which is clean and tight. A sensitive method 
of recording the radiation beam attenuation due to 
the oxygen absorption is also necessary. 


R. A. 


ment of a 


Evans' did the pioneer work on the develop- 


system which meets the aforementioned 


requirements. He developed a xenon gas discharge 
lamp (xenon has a resonance line at 1,470 A.), a vacuum 
spectrograph-monochromator to isolate the 1,470 A. 
wave-length radiation, and a camera for use with the 
ultraviolet radiation. He obtained some pictures of 
shock waves formed by transverse cylinders in a small 
two-dimensional flow field. 

The present program was undertaken to investigate 
the feasibility of the method for use in a large wind 
tunnel and its capability in vielding quantitative density 


measurements. 


EXPERIMENTAL APPARATUS AND PROCEDURES 


The studies were made in the No. 4+ wind tunnel of 
the University of California Low Pressures Project. 
The tunnel is an open jet, nonreturn type, operated 
continuously by a five-stage steam ejector system.’ 
Interchangeable axially symmetric nozzles are used to 
produce flows over a range of Mach Numbers. <A 
Mach +4 nozzle* was used for the tests described here. 
The tunnel facility includes pressure manometers; 
thermometers; flowrators; thermocouples with record- 
ers, necessary to measure flow properties; and a mani- 
fold equipped to handle bottled gas so that the tunnel 
can be operated on gases other than air. <A traverse 
mechanism inside the test chamber provides for re- 
motely controlled, counter monitored motion in three 
mutually perpendicular directions. 

The radiation source, spectrograph-monochromator, 
and camera described in reference 1 were modified for 
the present tunnel installation. The radiation source, 
an oil-cooled xenon discharge lamp with hot thoriated 
tungsten electrodes and a quartz envelope, was filled 
with a 90 per cent neon and 10 per cent xenon mixture 


ICAL 


SCIENCES FEBRUARY, 195% 
It was operated at 220 
he 


spectrograph-monochromator, supported rigidly in the 


to a pressure of 20 mm. Hg. 


volts a.c. with a nominal current of 1 amp. 


and 4, 


wind-tunnel test section (see Figs. 1, 5, con 
sisted of a vacuum chamber for the lamp with adjust 
able lamp support and a vacuum chamber with optical 
benches on which were mounted adjustable optical 
calcium 


element The optical system, of 


fluoride (CaF.), is shown schematically in Fig. 2. The 


supports. 


spectrograph-monochromator and the camera had a 


separate vacuuin system, shown schematically in 
Fig. 3. 

Five models were chosen for the flow pictures—a 
1°’ in. diameter sphere, a 90° cone-cylinder of 1-in 


diameter base, a l-in. diameter cylinder with axis 
) 


parallel to direction of flow; and two 72 
All models were 


wedges, 
one 1° 4, in. long and one 4 in. long. 
made of aluminum, and were supported by long stings 
attached to the tunnel traverse mechanism. All pho 
tographs with flow were made at a flow through th 
nozzle of 18 Ibs. Mach 
and a test chamber static pressure of ~114 u Hg 


hr. or at a Number of ~3.9 

A photographic emulsion, sensitive to 1,470 A. radi 
ation, which had good uniformity, high contrast, and 
low background density was required. Kodak Spec 
trum Analysis 41 U.V. Kodak 5S W R 
The exposure time for flow photo 


and (experi- 
mental) were chosen. 
graphs had to be such that for the photographic emul 
sion used, all portions of the negative were in the region 
of negative density, D, where the emulsion had a linear 
response. For the same flow configurations and emul- 
sion a much longer exposure time is necessary for pure 
oxygen than is required for air. the Kodak 
#1 plates were much less sensitive 


Since 
Spectrum Analysis 
but were more uniform and of higher contrast than the 
Kodak S W R plates, they were used for all the pic- 
tures taken with air flow. The Kodak S W R plates 
were used with pure oxygen flow where too long an 
exposure time would be necessary for Spectrum Analysis 


#1 plates. The S W R plates were exposed 90 sec 
The Spectrum Analysis #1 plates were exposed 50 to 
60 sec. Calibration exposures were made on each 


plate at appropriate pressures without flow through 


the nozzle. 
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VISUALIZATION 





Fic. 4. Speetrograph-Monochromator with covers off in tun- 
nel 1) Pressure tap lead tubes for jet-test chamber pressure 
balancing. (2 Model sting. (3 Nozzle. (4 Prism. (85 
Lens 43 mounting (6 Orifice #1 mounting rj Lens #2 
mounting. (8) Shutter 


The photographic negatives were studied by means 
of a Zeiss microphotometer modified to record on a 
Brown Function Plotter an indication of negative den- 
The 


microphotometer photocell output was fed into an 


sity as a function of position on the negative. 


electrometer and the amplified signal fed to the Brown. 
rhe recording roll was driven by a selsyn, with revers- 
ing switch, geared to the traversing drive of the micro 


phi tometer. 


ANALYTICAL PROCEDURES 


Mack? derived a theoretical expression relating the 
density of a gas in axially symmetrical flow (with zero 
density at the outer boundary of the flow) to the ab- 
sorption of a beam of radiation perpendicular to the 
axis of flow. In the Appendix a similar expression is 
derived which is appropriate fer an actual density meas- 
urement of a flow having a finite density, p;, at the 
boundary. From the Appendix, for a parallel beam of 


radiation of constant intensity, 


( r eR : 
i ae a[W(y), r] dy 1) 
PI Ted) 
where 
(p3 — pz), pi ; 
log [(100 — 1V3) (100 — IV, 
| Op: 


pi Ollog (100 — IVs 
div’ dy 
(100 — W)Vy? — r’ 


pis the gas density in the axially symmetric flow region 
ind is, therefore, a function of r, the distance from the 
axis of symmetry (see Fig. Al). & is the value of r 
at which the gas density becomes p;. p; is the value of 
the gas density in the uniform stream upstream of the 
flow disturbance. Subscripts 2 and 3 pertain to cali- 
bration values obtained under conditions of no-flow. 
L is the length of the radiation path in the test section, 


that is, the distance between the outside lenses #5 and 


OF 


L 


OW-DENSITY FL 


OWwWs QD 


#. yis the distance from the horizontal plane through 
the axis of symmetry (see Fig. Al IV’ is the reading 
on the Brown recorder in arbitrary units of divisions of 
Il” is, 
of the radiation intensity striking the 
Eq. (A-18 


the recorder chart paper. therefore, a function 
lsee Eq. (A-14 
photographic plate. £ is a function [sec 
of film characteristics. 

/: is determined by means of two or more calibration 
the 


density ratios are known and the corresponding read 


exposures where pressure and temperature o1 
ings Il and JV’; are taken with the microdensitometet 
IV and dIV dy are determined for a series of values of y 
For 


r, 8 is plotted as a function of y and the integration is 


from the trace on the Brown recorder. a chosen 


performed graphically. 

The analysis of the Appendix is possible only if the 
traverse, which is made perpendicular to the body 
axis, can extend into the region of known flow condi 
tions. The flow picture must therefore include a re 
gion of known undisturbed flow in the plane perpen 
dicular to the axis of symmetry through the region of 
unknown flow properties to be analyzed. Thus, there 
are limitations on the geometries of flow regions which 
can be analyzed with a radiation beam of given cross 
sectional area. 


EXPERIMENTAL RESULTS 


Some of the shock-wave photographs obtained aré 


shown in Fig. 5. Figs. 5(a), (b), and (c) are photo 
graphs taken in pure oxygen of shock waves formed 11 
front of a cone-cylinder. Most of the local variation 
in photograph density is due to the photographic 
process rather than to the flow configuration, although 
the body was moved to a somewhat different positio. 
in each case so that a different part of the shock wave 
could be seen. The effect of radiation path length is 
particularly noticeable in this case where the shock 
wave change in density becomes more apparent as the 
effective path length increases with distance from the 
cone vertex. Fig. 5(d) was taken in air, and the result 
of a relatively long effective radiation path (~ 1° , in. 
is readily seen in the high contrast of the region before 


the shock with that in back of the shock. Fig. 5(e) is 
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a similar photograph but of an axially symmetrical 
configuration. 
Fig. 6 shows a typical shock-wave trace with three 


calibration traces. The shock-wave trace was taken 


perpendicular to the axis of symmetry, at a location 


approximately 1/4 of the distance from cone vertex 


to shoulder. It was obtained from the same negative 


as Fig. 5(a). Ideally, for a uniform radiation beam, 
the calibration traces should be perfectly straight lines. 
Similar slopes at different absolute values would 


simply indicate a nonuniformity in the radiation beam. 
However, the 
shown is difficult to explain with the limited quantity 


somewhat strange emulsion response 


of data taken. Unfortunately, this response was not 
detectable either by the naked eye or a Weston Model 
#S77 densitometer and, therefore, was not anticipated. 
The limited emulsion response data made the inaccu 


racy too large for any density analysis to be significant. 


DISCUSSION OF LIMITING CONDITIONS ON THE METHOD 


In order to examine the effects of the physical 
parameters more easily, a simple two dimensional con 
Consider Eq. (A-2) of 
L), for 


The difference 


figuration will be considered. 
the Appendix for p = p,; = constant, (2% = 
each of two values of density p» and pz. 


may be written 


F i 4 = 
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p)CL = 
[(p 


In (J. [3) = (p3 — 


p2) — 1] (CL po R.T: ? 
where &, is the gas constant and /:, /; are the radiation 
intensities after passing a distance L through oxygen 
of densities po, p3, respectively. Cp is 350 em 

(see reference 5) for pure oxygen at N T P and 1,470 A 
C/R, 1s ~125 [°K./(cm.-mm 


radiation. Therefore, 


Hg 


and 


In (I / 125] Pz pe — | 
Fig. 7 isa plot of this equation. 

It can be seen that there are practical limitations in 
the detection of For a limited 
path length, if the density is low, In (J, /;) may be to 


close to the noise level of the system for good results 


resultant radiation. 


If the density is high, the signal may be too small t 
detect. 
if the photographic technique is employed. 


This would mean a prohibitive exposure time 
In the 
latter case, however, a gas mixture of a small quantity 
of oxygen and a large quantity of a gas with a low ab 
sorption coefficient, might be used; so that even though 
the overall gas density is high, the density effecting the 
radiation attenuation, the density of the oxygen, re 
mains low. 


For the S W R emulsion of Figs. 5(a), (b), and (c 


a negative density difference D, D; =~ 0.2 corre- 
= = aod { 
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Microphotometer traces of shock wave and calibration traverses, ~5.9 divisions per mm. of negative, scanning spot 150u X 
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| sponds to In (J; /3) ~ 0.5. This is a minimum level se T 
I) (2 for a good measurement with the S W R emulsion. p/p] 5 “ 
radiation This means that, for conditions described by the area 3 
1 oxygen | below In (2 J3) ~ 0.5 in Fig. 7, the present photo- 
(0 cm.-! § graphic process cannot be used. re % 
fe) + 4 | 
1,470 A 
Cim1.-1nm CONCLUSIONS 
Density measurements are possible by the ultraviolet 
radiation absorption method. The main limitations 
| on the method may be summarized as follows: (1 a ** 1 
Long effective radiation paths are necessary through ol 
ations in density fields of small gradients. (2) With a radi- " 
limited | ation source of finite intensity, the exposure time neces- kK 2 
y be too sary in high density regions with long radiation paths s 
| results is prohibitive. (3) The photographic detection of = } | 
small to} the transmitted radiation with present day commer- 
ure time | ¢jal ultraviolet sensitive emulsions requires that a large 
In the | number of calibration and flow photographs be taken to 
juantity : provide sufficient data for a statistical analysis. The P 
low ab- | statistical analysis is necessary because of variations 
1 though F jn emulsion response. This may not be a limitation 
‘ting the | when better emulsions are developed. _ 1s 
gen, re In order to achieve a system possessing the simplicity 
and reliability necessary in a practical wind-tunnel o 
and (c), instrument, modifications of the apparatus employed C 0.5 O 5 F 
2 corre- — in the present tests are essential. It is recommended Ba ssc ee: death 
that, before other studies are made, the following be = 

' investigated: (a) Use of a phototube instead of a Fic. 7. Flow parameters vs. radiation intensity ratio for oxygen 

camera for detection purposes. With a phototube 
"=+ | detector a lamp of greater intensity may be required. where ¢ = \V k* — y* and where p; is the known density 
Ses | 3 (b) Development of a fluorescent screen to facilitate in the region outside (yr > &) the axially symmetric re 

i optical alignment and focusing. (c) Construction gion. Of course, this assumes the radiation beam is 
T= =f | of remotely controlled positive positioning mountings parallel to the x axis and that the refraction is neg 
eames | {or the lenses, orifices, and prism, with verniers to per- ligible. 

; mit optical adjustment under vacuum conditions. Since 7° = x° + y’, it follows that for vy constant 
es = ' (d) Design of a simplified optical system having rdr=xdx = Vr? — y?dx. Whenx = 0,r = y, and 
Piast fewer optical elements. (The optical system con- when x = 3,r = R, so 
S==3 | structed included a spectrograph test section for lamp i ellen ah 

‘}-} | development purposes which is not necessary for the F(y) = -2¢ | at — Cp(L — 20 R? — y?) 
———+ | present application of the equipment.) ; rr = $F 

i A-3) 

\ | APPENDIX a 
=e Determination of Gas Density in Axially Symmetric Flow 4’ 

} by the Absorption of a Beam of Radiation 
ee r } 
tf | If a beam of radiation is directed in the positive x 
“+. | § direction (see Fig. Al), the radiation intensity J(x, y) 
= } | decays according to the well known law which may be 2... ' “oe iy) 

; expressed as . " 

; | ol,ox = —Copl (A-1) ‘ 
ae a ae . 7 

; i | where C isa known positive constant, p the gas density, 
+ | and x the horizontal coordinate. When the gas flow is 
a reas a . . . : » : = —_> 
t E | axisymmetric, the gas density is a function of r, the 
; =e radial distance. If radiation of intensity /o(y) is inci- 
7+ | dent on the region considered, parallel to x axis, then 
+ | /(y) will emerge and, from Fig. Al, 
InI(y) — In p(y) = -2¢| pdx — i . 
150u X J 0 


. — IF) = Fr) )) . . . . 
Cpi(L 22) a= f \y) (A-2) Fic. Al. Coordinate system for axially symmetric flow solution. 











A second change of variable may be made by putting 


S = R? — r? and y(S) = p(r). Then dS = —2r dr, 
S = R? —y? atr = y,and S = Oatr =R. Then 
0 y(S) dS 
F(y ( — 
J R?-y? V R? — y? — S 
Coi(L — 2 VR? - (A-4 
which may be simplified by substituting / Ke \ 
Then 
f(y F(y) + Co(L — 2V R? — y’) 
*t v(S)/t dS 
-c{ Gt A-5 
0 


‘ Vi-S 


Eq. (A-5) is of the Abel type® with unknown y(.S) and 


has the solution 


; 1 d 5% Git) dt 
y(S) = - aes (A-6) 
7 ( dS, o Rs —f{ 
or, after integration by parts, 
: 1 G(O) l *S G'(t) dt 
aU 4D tr Jo WS / 
where G'(t) = (d dt)G(t) G(O) f((R) = 0 
Now G'(t) dt = f’(y) dv and v R when ¢ = 0; y 
V R? — S? = rwhent = Sor 
1 ¢R f(y) dy 
pr 7 AS 
3 Oe Vy—r 
From Eq. (A-5) for Jj) a constant, 
l °R d d\ 
p(r) = . — in J -+ 
we dy Vy) 
| . R Coy dy ae 
3: ae Vy-rPVR-y 
After integration of the second term in Eq. (A-), 
l *R d dy 
p\’) -— pi = 4 In J (A-10 
7rC J, dy V/ yy? — 


Intensity /(y) was recorded by means of a photo- 
graphic emulsion whose density after development 
was measured by means of a microphotometer and 
Brown function plotter. The density of the negative, 
D, is detined by D = log(/; /;) where J; is the intensity 
of light impinging on the developed negative and / 
is intensity of light transmitted through the developed 
negative. For an emulsion with linear response D = 
I log J + C, = log J; — log J; (for a given exposure 
time). For linear plotter response IW = —A,/; + 
K» where VW is the signal from the microphotometer 
as read on the Brown function plotter, and A, and 
K» are constants. At infinite D, IV was set at 100, so 


100 — W = Ay, or 


log (100 — W) = log J; + log A, = 


log I; —T log I —_ C, + log ky, (A-11) 
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then 


(d. dy) InJ = [(In 10) T] [(dIV dy) (100 — IV 


A-12 
Substituting Eq. (A-12) into Eq. (A-10 
In 10 [ dW dy 
p\r - 1 = — ay 
wCr . (100 W)Vy-—r 
\-13 
From Eq. (A-11 
log (100 — IV) = —T log J constant A-14 


since in analyzing the negatives the incident light / 
was kept constant. But from Eq. (A-2), for p con 


stant, 


In J — InJp -CLp A-15 


Then, with In J) also absorbed in the constant term, 


log (100 Ws (LCT In 10) ps constant A-16 


log (100 IV’; (LCI In 10) ps; constant A-17 


Subtracting Eq. (A-16) from Eq. 


posing 


p p 
8 


(In 10) Cr : 
log [(100 — JI 100 I] 


j 


ap 


d {log 100 I> 


Now if both sides of Eq. (A-13) are divided by p 


l x 
p1 log [(100 i" 100 ie 
Lf div dy ; 
| dy (A-19 
™J, (100 — W) Vy e 


For /) not equal to a constant Eq. (A-10) becomes 


l ie d ad 
p(r pi — (nd In / 
me J, dy Vi 
A-20 
and Eq. (A-19) above becomes 
plr [(p pz) pil 
l = : ie 
p log [(100 — I] 100 — IV2)] 
L [ | dW dy dW.s/dy | dy 191 
©. 100 — W 100 — ir. | Vt gt ~ e 
where (dIIl2 dy) (100 — JVs) is the derivative of the 


logarithm of 100 — IT’ for p equal a constant. 
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An Introduction to Second-Order 
Wing Theory 
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SUMMAR\ 


In this paper the particular integral for the second-ordet 
iteration equation of inviscid supersonic flow, discovered by Fell 
ind the author,! is applied to wing theory This enables us to 
‘alculate second-order derivatives for wings of finite aspect ratio, 
thence to determine how far the second-order strip theory of 
Remarking that the most 


ind 
Busemann is in error for such wings 
interesting effect in second-order wing theory is the influence of 
thickness on lifting derivatives, and that corrections to strip 
theory will be largest on a wing with subsonic edges, we choose 
is our example the lift coefficient of a rectangular wing with a 


simple thickness distribution 


I\ 


Some details are given of the calculation, which is fairly com 


plex; use of the ¢ method introduced by Perry and the author? 
is found to reduce the labor involved very considerably 


INTRODUCTION 


i ie SECOND-ORDER strip theory of Busemann has 
been known for many years. If we apply this 
second-order theory either to an infinitesimally thin 
wing at incidence, or to a wing of finite thickness at 
zero incidence, we get a small contribution to the drag, 
but none to the lifting derivatives. To get a second 
order contribution to these derivatives, we have to 
consider a wing of finite thickness at incidence. Second- 
order lift is found to be zero, but, for wings of practical 
thickness, the second-order correction to the first 
order pitching moment may be very large. 

It would be of great interest to know whether the 
predictions of second-order strip theory are any guide 
to the behavior of wings of low aspect ratio, but this 
requires a full solution of the second-order problem, 
which has not previously been available. Martin and 
Gerber‘ have calculated the damping in roll of an in 
finite sweptback wing, and Clarkson’ has made a pre- 
liminary study of the general problem, but no second 
order derivatives for wings of finite aspect ratio have 
been published hitherto. 

The rectangular wing has been chosen as an example, 
since on a wing of simple planform, the ‘‘aspect-ratio”’ 
corrections to the second-order strip theory lift will 
probably be small (the corresponding first-order cor- 
rection is zero). The thickness distribution chosen 
seems to be the simplest which gives zero thickness on 
all the edges (the theory would not be applicable to 


wings with thick edges). 
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FORMULATION OF THE SECOND-ORDER PROBLEM 


The systematic development of second-order theory 
is due to Van Dyke,* and this section should be read as 
a sequel to his paper. Our notation is similar to his, 
but differs from it at one or two points. Equations 
given in Van Dyke's paper are quoted here without 
proot. 

Ward’ has shown that, in first-order theory, the prob 
lem of a lifting wing of finite thickness can be sepa 
rated into an antisymmetric or lifting problem in which 
the thickness of the wing may be ignored, together with 
a symmetric or drag problem in which the wing may be 
treated as being at zero incidence Let denote 
the height of the wing surface above the midplane, 
while Z(xv, v) denotes the distance between the mid 


plane, and s = O, which is taken as the “‘plane of the 
wing’ (the supersonic free stream being in the x dire¢ 
tion). Then 

g OY On 


is the thickness incidence, while 


a OZ On 2? 
is the lifting incidence. The overall height of the wing 
surface above the plane s = 0 is 

S y+Z 


depending upon whether we are on the upper or lowe 
surface of the wing. Let the potential of the anti 
symmetric problem be @ and that of the symmetric 


problem be @. The overall linearized potential will be 


sy 6+ 3 
on the lower surface, and 
®,, = @ @ }) 


on the upper surface, the functions 6, @ being defined 
initially in s < 0, and being extended to s > O by their 
symmetry properties. The second-order potential will 


be of the general form 
r= (6 F+ @ 


(this is not an equation for ¢ The antisymmetric (or 
lifting) portion of ¢ is of the form 6 X 9, the terms 
6°, ¢° merely contributing to the drag. 

It is convenient to evaluate the potentials and their 
derivatives on s = 0 rather than on the wing surface. 


If we follow this convention, the second-order pressure 
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Fic 1] rhe thickness distribution and dimensions of the wing 


studied in this paper 


coefficient is given by 
C,’? = —2¢, — 22,%,. — ®,? — %," + B°® ») 


where ® denotes ®, or ®, and ¢ is the second-order 


potential, which is given by 
g = I’, + Y+ oa 6) 


Here y, the “irreducible part” of ¢, is a particular 


integral of 
pg = Gy a Vaz 


and the complementary function @ is a solution of the 


~ BY. = 2M%82N@,6,, (7) 


linearized equation. @ is, of course, determined by 


the boundary conditions, which are 


@, = (0/Oy) (2,%,) — B°(0/On) (2,,%,) — 
A[*°@(0°s,. Ox?) — 


on the wing (Sa) 
@= 0 off the wing, (Sb) 


As in the corresponding first-order problem, Eq. (Sb) 
is a consequence of the antisymmetry of the potentials. 
It follows from Eq. (S) that as long as the subsonic 


edges of the wing are independent —1.e., as long as 


A’ = 28b/c > 1 


@ is given by Evvard's theorem (see reference 7); the 
dimensions 0, c of the wing are defined on Fig. 1. 

The explicit calculation of a(x, y) from Evvard’s 
theorem would be a formidable matter, and, fortu- 
nately, it is not necessary. The quantities we actually 
want to calculate are integrated derivatives of the form 


Q= 2 | | f(x, y)C, (x, y) dxdy 


wing 


where the distribution function f is 1 for lift, x for lead- 
ing-edge pitching moment, and so on. We shall con- 
fine ourselves to the case of lift (f = 1), but the general- 
ization to more complex forms of f will be obvious. 
Thus in dealing with the complementary function, we 
wish to evaluate integrals of the type 


ff @,(x, y) dxdy (9) 


wing 


Now @ satisfies the first-order equation and, therefore, 
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Q(X, Vb," (x, v) dxdy if aS Ya 


where ¢‘*’ is the linearized potential induced by a re 

versed flow over the wing at constant unit incidence 

Since @, is given on the wing by Eq. (Sa), Eq. (9a 

enables us to evaluate integrals of the form (9) without 

recourse to Evvard’s theorem. In 
We now substitute Eqs. (3) and (4) into Eqs. (5), tec 


6), and (Sa) and separate out the antisymmetric part 
of the second-order pressure coeflicient: we als 
change our notation slightly, introducing reference 


incidences ay, &, and writing 





Eq 
a(old) = ay X a(new | 
' SCC 
and similarly for €; the new functions a, & are non- pr 
dimensional functions of order |. 
The quantity we are actually going to calculate is 
> > ill 
by at 0°C, Oa Of = rs My c ae y dxdy 10 
where C, 1s the lift coeflicient, and 
r= 1/aoé-S 
S being the wing area. Using the first-order boundary dit 
ie lon 
conditions, 
in 
U & £ @ r 
we may show that b A 
Cua * + Lb ZS 
In this equation the 
l ae id are 
“2, the 
stit 
where 
h, = —4aé, h. = —4) (0a On), 
all 
h; = —4Z(O€ Ox) (12a) | 
7 4 P : 7 
I, =! I(x, vy) dxdy (13 
wh 
where 
1, = +467¢,0,, ly = —46,0,, 
lz3 = —412°(0 Ox) (6,6 
(13a In 
iy = —4.11°(0 Ox) (6¢,), lb = —dyr, the 
@= LY, free 
and 
0.1.4 
J,=T ff o,"®)(x, y)7-(x, vy) dxdy (14 
wh 


wine 


l Qa 
by a re 
cidencs 
Vq. Ya 
without 
“qs 5 ; 
ric part 
ve als 


elerence 


undary 





} 
}” 
l2a) | 
(15 
(13a) 
(14) 
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st 


where 


bae(O Ov ZO 
tO~-a@ O Ox ZO 
1 l/*ap-6(Oa@ On 


i, is = AY 


In Eqs. (13a) and (I4a), Y and y are particular in 


tegrals of the equations 


O*y, = 2.1/°8°N6,,o lja 
OY 2.1/°3°.N6,.0 15b 
Eqs. (11) to 15) constitute a closed formulation of the 


second-order antisymmetric problem. In our present 
problem, we have 


a l, S Ihe 


and therefore, since } is zero on all the edges, 


H, = fH. = 0 HT; = 4 


DESCRIPTION OF THE ¢ METHOD 


The terms J, and /, can, in principle, be calculated 
directly, but such a calculation would be extremely 
long and tedious. The ¢ method, which was introduced 
in references 2 and 3, simplifies this work very greatly 


he first-order drag potential is given by 


(x, vy, OU) = 


; . 5 ah si 3 
& re ax dy 
E(x", v 
7 . Vix — x’)? — Biv — vy’ 


the integration being over the forecone of (x, v we 
are not interested in the values of the potentials off 
the plane s = 0. In this equation, we make the ¢ sub- 


stitution 


y=y— f(x — 2x’ 16 


+ el p 
so = (ls = 
KS, 3; 0) = = EC eee A Ly 
1/3 | J 


where 
Px 
2i((; x, y) = E[x’, vy — ¢(x — x’)] dx’ IS 


In this definition, £ is to be understood to be zero off 
the wing. These expressions may be differentiated 


freely under the integral sign—for example, 


& 71/3 : d¢ 
0,(Xx, y, 0)= — &(C; X ¥) 
a « 1/3 


where 


ws 
l| 

It] 

- 


Os 
. Ox’ 


From Fig. | we see that on the wing 


£ 1 — H(mnx y b H(mx + 4 h INSb 
m being the slope of the bend lines AC, BD (.1 
3m the origin is taken at the midpoint of the lead 
ing edge.* Substituting this into Eq. (1S), we find 
= v(b + vy — fn 
mx + y b) (m + €)| H(ma 
mx — ¥V b) (m C)) ) (my 
He) 1 by-F >O0 (19 


a formula which will be used repeatedly in subsequent 
work 
The case of the lifting potential is slightly more com 


plicated. The general formula 


is still valid, and ¢, is still equal to aa on the wing 
however, outboard of the wing, but inside the Mach 
cones of the points A, B, ¢, has a non-zero value which 1s 
initially unknown. This value can be determined by 
conical flow methods, for example, and if we then 


define 
1,(¢: x. y) = o.1. c(A \ ax 20 
we shall have 
l iis a4 
g(x, Vv) = - Pilge a 2 
we “ Vv i +~4 


The matter is discussed at greater length by Leslie, 


who shows that 


A,=x-— jx b+y) ¢]) X 
\ } \ /p 
Myx — (6+ y) ])} Va + Be) Be (20a 
for ¢ > 0; there is a similar formula for ¢ < 0 


CALCULATION OF /, TO J, 


The calculation of these items is quite straightfor 


ward; we shall take /; as an example, since some ot 
the work will be useful lateron. Now 
.» 
I, = —4)/°T 6.(Cc, viglc, ¥v dy 
and we therefore define 
iia f ; - 
O3(C, K) 41} oe fi(K: ¢ Vv) ay ZZ 
; Phe « b it ; 


* Here (x) is the unit function, defined by 


H(x) = 0 if «vo 
H(x = | if x 








oo 
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/;,is then given by 


1/p 
. 


1M? {| ’ Bde Bdk 
= — — 53(€, «) 
TP J V1 — Bo? V1 — Br? 


1/p s 


Since a and £ are both symmetric functions of y, we 


have 


—x) = S3(¢, x S3(f, — «) = S3(—¢, x 


S3(—¢, 


from which it follows that 


Sir [ 8 Bdk 
, 


rp? JO VI — Bx? 


adel ss Bd¢ 
93(f, K) (23a) 
1/6 V1 — B¢t° 


We therefore only need to know =, for « > 0 [ef. Eq. 


(19) }. 
Regarded as a function of ¢ for fixed positive x, S: 
changes analytical form at 


¢c( = —(m — kK c= J Cc = kK 


If 1’ < 1, there are further changes of form at ¢ 

tm, but our theory depends on Evvard’'s theorem, and 
is therefore only valid for .1’ > 1. The change of 
form at ¢ = —(m — x) occurs when the line of ¢ dis- 


continuity 


2, and /; must therefore 


A" = 2. If we 


This can only occur for .1’ 
have a change of analytical form at 
confine ourselves to .1’ > 2, as we shall do in the pres- 
ent paper, S; has differing analytical forms in the three 
regions 

—(1/8) <¢ <0 Os € <x c< ¢ < (1/5) 


The calculation of these three forms is quite straight- 


forward: substituting them into Eq. (23a), we find 


B21;/M? = 4[(2A'/V A’ — 1) — 1] X 
{1 — (1/2A’)] — (4/77A’) (4K. — Ki + Ly + 2L3) 
(24) 
where 
K, =K,(A’) = 
iis u” : ) du : 
. at V2 (25a) 
Jo Ao — u Vi- wv 
and 
(25b) 


F oT. [i u” du 
= Se Boetd = p 7 
JO A -— uUVu(l + Un) 


These quantities satisfy the recurrence relations 
. nom 
K, =A eo 1 — Age 
Ri _ A'L,-1 ica Dy, 1 
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where 
| 
>] . 
: b au 
A, = “usin Vu = 
70 Vil—-u 
n— | 
\ D, + D 
n ”) 
and 
o u" du ae 2n 
D, = . D 
“YO Vull + iu a 2H 
The A, and L, have been computed by using the re- 


currence relations, the ‘“‘down-ladder’’ process being 
started at values of » high enough for the asymptotic 
formulas 


7 7 . , 
; Sin 4 dé 
2(4° — 1 0 


to be sufficiently accurate (for example, a 25 per cent 
15) only affects Ls (15) to one unit in the 
I 


It 


change in Ls 
fifth place). The auxiliary quantities \,, D, can 
obtained by an “‘up-ladder”’ process from 


1.53196, 


Do = 21In (1 + V 2), \ 


(A) has been calculated by numerical integration 


Values of A, and L,, for suitable ranges of » and 4 





have been computed by the above process; with these 
values, we can calculate 6°/; .\/* as a function of Al’ 
It should be noted that Eq. (24) is only valid for .1’ > 
2. Tables have been prepared giving the values of 
K; Log Xn, and D,, S-(¢,.«), Rx) (these quantities 


are introduced later) and the /, and //,. 
in the JOURNAL has prevented the inclusion of these 


Lack of space 


tables, but the author has a few copies available for | 
distribution. 

I;, Iz, and J; can be calculated in exactly the same 
way. Wedefine 


> 
ANG; X, WE\K; X, V 
be J e : 


Win 


dxdy 


| Ot, .Of1, _ 
ti s9 K; xX, v) dxdy 
Pho J J Ov Ov 


in 


So(f, kL); = 


and 


. 
I c ; 
aylC; C, V)S1(K; Cc, y) ay 
} 


S4(C, Kp= 
Pbc e 


where 
a, = 0.1;/Ox 
We then have 
Sf, [ Bdk B dg 


+] 
7 | Bs ae 
0 V1— Bx. 1/B V1 — B7¢? 


rT: 


whe 
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} wher where 
f, = —(1 8 fy; = fy= —(M? 8B SIZN 3°x7- 3 dh 
= x 
m3? Jo V1 — B% 
CALCULATION OF I; AND I, . F 
, 4 
~ D S CK Og 
In this group, we shall take as our example ‘ 1/1 — gt 
, SIPN f' 
I= —40 | ¥i(c, y) dy / ox Vv 2 dx & 
: r*B> JO 
D - as ii 
It follows from the results of Fell and Leslie! that 23\5 9 * adel se Ob 
: . ¢ h \ ie 4 
! 7 De? , y 
the re- y » vy ~ = A N 7 SISS (Oa ON) x 
s being v1, V1) (O°E Ox: v2, ye) X In this last equation 
mptotic O(X, V, 35 X1, Vij Xe, Vo) dxydxedyid 26 g3(f, K) = S3(f, «) + K(O OK) S3(C, & 
where Evaluating the integrals of Eq. (29), we find that /; can 
be expressed in terms of the A, and L,. Js may be 
; _— R-—R):(R-—R, -—- R-R R R treated in the same way; its relationship to S; is the 
i -" (R-—R,)-(R-—R.) +/R—R,|R—-—R same as that of /; to. S 
13-X ba 
CALCULATION OF THE TERMS ./; TO J; 
er cent \. being the quantity introduced in reference 1. We _ ; 
t in the substitute into Eq. (26) the standard ¢ relations Phis group of terms oS ee straightforward, and as 
can te an example we shall consider 
= (1s — Hy and Yo = yv—xK(x — A se 
: : : 4 O 
' : J tla Q vr, V Z6 axd\ 
ind that equation becomes on = = 0 (we are not inter ve O' 
ested in the value of y, elsewhere a 
In reference 3, itis shown that 
eer pee in ££ l dé 
ation ; [?, oti. ak delle 30 o _ Vy 30 
ind 1 Dn? J. T « V1 32 
1 these , 
a : where 
ol F Oa 
5 ‘i > \ = 2 is \ - (1s = Vy dn x Vy ne as ¥ l if D + ; l + OC x 
f ‘ On yal 
ues ol aa. =Pj VSO 
vate ; O°é s 30 
ntities x—-xXx a ae ee dn IT »{(1 O¢ ( \ 28; . 
space . On» 
these ; 0 otherwise 
- | Where, from Eq. (26a 
le for | - ' 
ole Chen if we define 
1 — B%x) — V1 — By? V1 — Bx i 2 
- Same 7(¢, k) = log Rif, k Vict 4 x 
(1 — B-tK +Vi1—-— Ber VI — °K" Pbc J « OV 
We may readily show that Eq. (27) 1s equivalent to I, = { ‘ 
- 4 \ y AK \ \ ad vd \’ ol 
7 l/s Ov 
IN (¢ O° we ; 
Wilx, yy = oc. & “ We have 
2a" 6s OCOK 21/3 
1/3 S a) dk 
CRAVE; XY, WE KS XY, Vv) | d&dk me VI — ~*« 
ind, therefore, that Rie d¢ ») 
1/3 = “1? V 1 3-¢ 
2M2N ( O 
] es a(C, k R, has the same symmetry properties as S;). In our 
- he OF OK particular case, Z = x and using Eq. (30a) we have 
[cKS3(F, «)] dfdx (28 i ff wh... 
R, = : x axay 
. a = 4 “ ; POC J eo Ov 
S; being the quantity defined in Eq. (22). Integrating DMNC , 
Eq. (28) by parts, and using the symmetry properties | Z i O=1 . 
yy Z y proj . FP 
—— a - ZA 
c Ol Ss, we find Phe \J Ne J MD ov 


i Ip = Isa + Ip 29 the lines MD, NC are shown on Fig. 2). These in 
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1(+BC x i BC)(e_x 
™ s c ae | ld co - a - 
- 2p - Be 
A M N B 
Paneer eet = hee 
/ ‘ae | 
Ww | 
x | 
} 
| 
} 
| 
| _—_—n J 
Fic. 2. Diagram illustrating the area of integration for R 
tegrals may be evaluated without difficulty. [From 


Eq. (32), it will be seen that a change in the sign of 
¢ in any of the terms comprising X, leaves /; unaltered: 
this property can be used to simplify the formulas 
for Rk; and the other /#, defined below.| Substituting 


the resulting formula into Eq. (32), we find 


k= (GA + 1DWV2)/(A’G ' + 1)?] 


a formula which is valid for A’ > 1. 
The terms /, to /s can be dealt with in the same way. 
We define 


| 2 
R,(¢, x) = {| Vic; x, y)re(w; x, vy) dxd} 
; De « « 
where 
ry = (O Ov) [¥ (0-1, Oy) |, rz = (00x) (Zé), 


rs = (O Ox) | Vas), 
rs = =,(Oa On), re = -1,(OF On) 
The required terms are then given by 


8d¢ 


Ng, »l 2p B dx 21 8 
RAg, «) 
J0 1p Vi-— i 


V1l— BK" 


ot 


where 


CALCULATION OF J; AND Js 


In this group, we shall take as our example 
J; = 4r o,\"' (x, v)Yi2(x, ¥, 0) dxdy 
¥i- can be obtained by differentiating Eq. (26) under 


the integral sign, to give 


viz = (MPN 22?) SfSF (Oa Ox) (m1, y1) X 
(O°E/Ox2”) (x2, V2) (O0,/O8) (x, ¥, O; X1, Vi; X2, V2) X 
dxdxedydy. (33) 


Differentiating Eq. (26a), we find 


Oa 27 Vy ~. 


Oz x — x) (x - 


LV = Vi y~- vy Sy — Vy 
él - ot 4 : Bas 3 5 
= Ls x= XL G = Ds 


1e Dirac delta function 6(.) being the derivative of the 


— 


t 
unit function //(x). We now substitute from Eq. (33a 


into Eq. (53), and, in the resulting equation, we set 


Ihe delta factor of Eq. (33a) then becomes 46(x 
which shows very clearly how suitable the ¢ method 
is for this type of work. Carrying out the y integra 


tion, we have, finally, 


IN (' 
V1- c, % Lh Y, V V 1 I-K- Th 
where 
i [ Oa 
glx: xX, ¥) = VY) — NX “1, V— K(X — Xx x 
: Js On. ’ 
Org 
\ \ is = 2 AN AN ) 
ON" 
We then define 
K;(¢, x) = Mics x, vyulk; x, v) dxdy (36a 
Pbc Je ‘ . 
and J; is then given by 
SIZN ('? 
J; = SNE BWV 1 — 8°" dk x 
r-6- J0 
71 5 , F 
Oo d¢ 
R-(¢, « 36b 
v-ls V1 — 8° 


The evaluation of u requires a certain amount of care 
Substituting from Eq. (1Sb) into Eq. (35), we find 


Oa 
“= HT'y +b — xx Xx; = Xs 
J. Ox; 
(x1, ~~ Ce = & 0 (Xx dx\dx>» se 


m ‘ ‘ ei Oa 
( ) H(y + 6 — xx Xv}; — X x 
m— kK J .« ON) 
is 7 b - KX . 
[x1,¥ — K(x — x) /6’ (xe — dx\dx. — 
we — « 
( m ) [ i Oa 
XY) — Xs 
ow « Ox) 


We = « 


Kv —y+b 
ae : rire 
[x1, 9 — K(x — xX) 16’ [xe — dx \dx» 
we «K 


In the region of y + 6 > xx, we have 


a= F(x, ) 
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so that the first term is equal to “05 Med 
son 0) 0 7a 
nd the second to 
a i Hiy + 6 — xx) |}Himx ; 
F -O3 
| V6 | (7974 v — 5b) (m — x)}} 37b 
ive of the Ss 
. - ° - A: > 
Eq. (33a [The third term needs more tormal treatment we find aoa Me 
Vie a 
we set that its value is Mr25 
xX: } 7 & oly = 3 O Or i T - 
- ¢ HI \ +t INAX Qa vi 
kK E . 
method \dding together Eqs. (57a), (87b), and (37c¢), and sub : : 
. ; age “ee co ae 5 2 3. 4 € é K S 9 iO 
integra. | stituting into Eq. (56), we obtain formulas for R; and A 
thence /;, which is found to be zero; we have not found FIG Graph of CL, ag aga 1 f V2, 2 2 
ny basic explanation for this 
ae /. mav be treated 1n exactly the same Way; we find since J,+/ t 
j SMAN : : . Tables of /;, /2, and for a selection of values of + a 
} 3 3dx X . * 
} r3? Jo have been prepared by the author Using these tables, 
“1/4 se Cr. 4 , can be computed for any value of 4 and MV 
O d¢ ae —_ . 
R3(c, x In Fig. 3 we show graphs of ¢ »¢, aS a function ot 
: \ — pc . ; 9: 
: lL—, {for UM = V2, M = 2, and M = 2.5 (we take 4 
ie 1.4 
a From these graphs, we see that C; ,; is always 
! > . . . a? 
\ ) R, = - Vic: x, vivlk: x, v) dxd\ negative, and that it tends to zero as 1’ tends to in 
aii finity, in accordance with strip theory. The maximum 
ad absolute value recorded is 0.47 at .\/ V 2 and 4 = 2, 
but since the strip theory value of Cy, ,¢ is zero, it is 
364 < % ao rather difficult to say what this means in terms of de 
mor | mn eee ee eS x parture from strip theory. In this connection, w 
may note that, according to strip theory, the second 
Og — order lift of the panels AOB, COD (see Fig. 1) is 
Xe, V — K(X — Xo)] dxida 
Ox ACr ee = +(M2N — 2)/8? = +28 
By meth similar > -d for > analysi ; 
31 ethods similar to those used for the analvsis of at M = V2: the maximum absolute value recorded 
u, We mav show thx ‘ NE 
26) , ses above is about one-sixth of this. 
AUD - . 
-[m (m — «)|H(b + y xx) H(mx — y — b) - Rather larger values of C; .: may be expected in 
i ae La ae \ extenc ; = 
m'(m + «)) H(y + mx — b) {(8/O) (kay) — the range l = A’ < 2, and the extension of our results 
am care. | , ) to this range is now in hand. Further work is needed 
y+mx — 0) (m+kK Oa, ON), : ; ; 
id y only on the /, terms, the formulas given for the /, being 
from which formulas for A; and J; follow at once. valid for 4’ > 1: the main work involved is the cal 
Che formulas for /, obtained in this way are all valid culation of the S, for << — (m — k 
lor to? 1. The method of Bonney* can readily be applied to the 
’ . ° 
wing under study. It gives 
i DISCUSSION OF THE RESULTS 
BC, az = (APN — 2), [A/V AA’ + 1 
Combining Eqs. (12), (13), and (14), and substituting 1+ (1/24! 
1/2A‘')} 
' into Eq. 11), we see that the second-order lift coeffi- 
> cient may be written in the form According to this formula, Cz 4,¢, 1s small and positive, 
ia : ; its value at 4 = 2and VJ = V 2 being 0.186: Bon- 
°C, at = M*N1,(A’) + M*le(A’) + Bls(A 38 . WH, 
ney’s method is therefore not satisfactory for wings 
Where /;, /2, and /,; are functions of 1’ only, defined by with sloping bend lines. 
| h = (8 IPN) Js + Je + Jz + Js 39a ' 
CONCLUDING REMARKS 
lo = (8? AM?) (3 + 1h + Js + SJ 39b 
Although the strip theory prediction that second 
sa and order lift is zero is not correct at low aspect ratios, it 
l = 4 + I, = i —- Ji — Js = J, — J; 39¢ *Because of space considerations, these tables have been 


omitted from this publication. However, copies may be obtained 


=J+1,+/2+ Js from the author 
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does seem that the second-order correction to the first- 
order lift is unimportant for practical wing shapes. 
This remark does not apply to pitching moments, and 
it would therefore be of the greatest interest to know 
whether the strip theory value of the second-order 
pitching moment is seriously in error at low aspect 
ratios. The calculations reported above suggest that 
this may be so, but to check the point it is necessary 
to calculate the accurate pitching moment coefficient 
about some suitable axis. 

It will be noted that these methods could be applied 
to the calculation of accurate quasi-static damping 
derivatives; the corresponding strip theory is due to 
Van Dyke.’ 
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APPENDIX 


We shall now show that when the incidence is inde 
pendent of y, our Eqs. (26) and (53) lead to the stand 
We shall use a one-incidence 
In this 


ard strip theory result. 
one-potential formulation to simplify the work. 
formulation, the irreducible part of the potential is 
given by 


me 


da 
dx, 


d : a 
dx." 


ul 
dx\dx» | oA(X, ¥, 3; X1y My Xe, Vo) dyidys 
Integrating by parts, we have 
x Bz 
M?N da da 
= = f(x, ¥, 2; X17 Xe) dxjdx 
t}r- Je dx; dx» 
where 
Vach 
“ "( (0c, Oa 
=> = dy\dy. = 
i Je Ox; Oa 
C one 
T- 
= (2 — xX) — X) 
Thus we have 
y = —(\?7N 267) [Bzi a(x aad Bz) }? — 
a(x — Bz)Z(x — Bz)| (40) 


= (1/2)1?N4,[sb, — | 





which is identical with the formula given by Van Dyke® 
for this case. 


We can now go on to consider Eq. (35), which, for a 
one-potential two-dimensional problem, reads 


x 
w= ff Xi — Xe 


Substituting this into Eq. (34), we get 


, da 
dx,dx. = Z(x) 
dx 


d*a 


dx»” 


da 


dx) 


¥.(2 = 0) = (M?N/28)Z(da/dx) 


in agreement with Eq. (40). 
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Supersonic Flutter of a Cvlindrical Shell’ 


JOHN W. MILES* 


/ niverstly of California, hos Angeles 


(1) SUMMARY 


Phe aerodynamic instability of traveling waves on a cylindrical 


shell exposed to an external supersonic air stream and containing 


an internal fluid is examined on the hypothesis of short wave 
length (compared with the radius and length of the shell The 
stability criterion for an empty isotropic shell and the negative 
damping ratio when this criteria is not met are determined. It is 


shown that the negative damping ratio calculated in the absence 
of structural damping decreases monotonically with an increase in 
structural damping, but structural damping cannot prevent in 


stability and may render an otherwise stable motion weakly 


unstable 


(2) INTRODUCTION 


_ FOLLOWING ANALYSIS of the dynamic instability 
(panel fluiter) of thin-walled cylindrical shells 
under the action of an exterior supersonic flow is pre- 
sented as a sequel to an earlier formulation of the two- 
dimensional problem.' In that paper panel flutter was 
studied as a traveling wave phenomenon, following 
Kelvin's classical analysis of the surface waves pro- 
duced by wind blowing over water. We remark at the 
outset that all of the essential features of this approach 
to panel flutter are elucidated by the two-dimensional 
analysis without recourse to any approximations beyond 
those implicit in the assumptions of small disturbances, 
a perfect fluid, and an infinite surface; however, further 
analysis and additional approximations are necessary in 
order to render the results applicable to practical con- 
figurations. 

These principal conclusions of this two-dimensional 
analysis were as follows. 

Waves traveling upstream relative to the panel are 
necessarily stable (this is significant in any attempt to 
extrapolate the results to standing waves). 

The critical air speed, l’., above which flutter is pre- 
dicted, may be related to the speed of propagation of a 
given wave, |, the mass parameter, (p,A m), and the 
Mach Number, .J/, according to 


U. = Vol(A)f[(oar/m), AM] (2.1) 
where \ is the wave length, p, the air density, and m the 
areal mass density (/ was determined implicitly and 
presented graphically in reference 1). 

l’. exceeds J) for all values of the mass parameter 
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t Many of the details of this further analysis, especially as re- 
gards justification of the approximations, are omitted in the 


present paper but may be found in references 2-5 
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(p,\ m) and Mach Number .\/ but is asymptotic to |’ 
as either of these parameters approaches infinity. 
If pyri< 1, the instability approximates a 


simple (exponential 


mand J < 
divergence, and the critical speed 
is given by (assuming air to flow over only one side of 


the panel 


ai 2am r)V (2.2 


If p,A\ < mand .\/ > 1, the instability takes the form 
of a slowly divergent oscillation, traveling upstream at 
the velocity of sound a with respect to a point fixed in 
space, and the critical speed is given by (im two dimen 


S1OHS 
U.= Veta (2.3 


The foregoing criteria imply that planel flutter is 


critical in the (spectral) neighborhood of minimum 


wave speed (I). The mass parameter for cylindrical 
configurations usually would be sufficiently small in this 
neighborhood to render subsonic divergence, as gov- 
erned by the critical condition (2.2), unlikely, and we 
shall consider further only supersonic flutter at small 
values of the mass ratio. 

The most striking new aspect in the flutter analysis of 
a cylindrical shell, relative to that of a flat panel, is that, 
in virtue of the curvature of its surface, any deformation 
of the shell is necessarily accompanied by direct exten 
sional stresses for which the variation of wave speed with 
wave length is opposite to that associated with flexural 
stresses; the result is a true minimum in the wave speed 
analogous to that found for water waves in consequence 
of the opposing effects of gravitational and capillary 
forces. This minimum wave speed occurs at a wave 
length that is of the order of the geometric mean of the 
radius and thickness of the shell, thereby justifying the 
approximation of an infinitely long shell. The smallness 
of the wave length compared with the radius of the 
cylinder suggests further that the aerodynamic forces 
are approximately independent of the curvature of the 
surface and essentially similar to those associated with 
a traveling wave on a flat panel, and that they may be 
calculated as perturbations on the /oca/ steady flow 
i.e., the values of p,, a, and LU’ may be taken as those in 
the local flow rather than at infinity.’ 

A second, important difference between the two-di- 
mensional flat panel and the three-dimensional cylindri- 
cal shell is the necessity of prescribing two parameters, 
say wave length and helix angle, to specify the three- 
dimensional wave forms. It follows from this considera- 
tion that the critical requirement that the wave move 
upstream with sonic velocity relative to a fixed point in 
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Fic. 1. The surface wave ¢ and the relative air speed lL’ in the 


local coordinate system (x), 91, 21, ¢ of the shell 


space may be met by a (finite) continuum of waves in 
the two-parameter spectrum, whereas the condition (2.3) 
defines a unique wave form in the two-dimensional 
problem. We emphasize that this follows from the in- 
crease in number of dimensions, rather than the change 
in shape of the panel; in particular, a two-parameter 
spectrum of waves also exists for a flat panel in three 
dimensions. 

The effects of the fluid below the panel (and at rest 
relative thereto) were neglected in the two-dimensional 
analysis, but the fluid enclosed in a cylindrical shell, 
especially if under pressure, may assume a dominant 
role in respect to panel flutter. Pressurization leading 
to relatively large, static stresses will be important inso- 
far as the membrane stresses associated with variations 
about these static stresses (and proportional to the 
curvature induced by the traveling wave) become 
comparable with the extensional and or flexural stresses 
induced by the traveling wave. The importance of the 
inertial effects of the internal fluid will be proportional 
to the parameter p;A /m. This parameter usually would 
be sufficiently small for an internal gas to render its 
inertia negligible compared with that of the shell, while 
for a liquid it could be sufficiently large to render negligi- 
ble the inertia of the shell—in the latter instance the 
flutter problem may be regarded as an extension of 
Rayleigh’s analysis of the disintegration of a jet,® 7 with 
the elastic forces exerted by the shell being regarded as 
analogous to surface tension (the analogy being par- 
ticularly apt for those stresses associated with the pres- 
surization). 

The compressibility of an internal liquid would be 
negligible insofar as the flight speed (and, therefore, the 
critical wave speed) were small compared with the ve- 
locity of sound in the liquid; the latter criterion might 
not be met by a gas, of course, but then the effect of the 
internal fluid probably would be negligible in virtue of 
the small value of pA/m. The possibility that an 
internal gas might prove important if the frequency of 
panel flutter were to coincide with one of the resonant 
frequencies of the cylindrical cavity is ruled out for the 
more important cavity modes (for which the wave length 
is of the order of the radius) by the relatively small 
value of the wave length of minimum wave speed. 

Structural hysteresis was neglected in the analysis of 
reference 1 and was not considered explicitly in ref- 
erences 2-4, although it was erroneously stated in ref- 
erence 2 that small amounts of structural damping 
(expressed as per unit diminution of amplitude per 
cycle) could be directly superimposed on the negative 


damping ratios calculated in its absence. A more pre- 
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cise examination of the neighborhood of maximum in 
stability reveals that, while structural damping is bene 
ficial in this neighborhood, no finite amount thereof can 
prevent the instability predicted by /imearized theory 

It might be thought that, in analogy with the action 
of oil on disturbed water, the viscous friction between 
either the external or internal fluid and the shell could 
provide appreciable damping of panel flutter. We find, 
however, that the damping per cycle on the assumption 
of /aminar flow is negligible compared with that which 
might be expected from structural hysteresis. It is quite 
conceivable, on the other hand, that the amplitude oi 
flutter might be limited by turbulent damping. 

A more complete understanding of panel flutte: 
even in a linear approximation—demands an investiga- 
tion of the refractive (as opposed to viscous) effect of the 
steady flow boundary layer. We may anticipate that 
the existence of a region having a thickness not small 
compared with the wave length and a steady flow ve- 
locity appreciably less than the flight velocity could al- 
ter the results obtained herein in a decisive way. This 


question is now under investigation. 


(3) THE EQUATIONS OF MOTION 


We require the equations governing the propagation 
of a transverse wave over an infinite shell that contains 
fluid internally and moves with velocity L’ relative to 
the air external to it. We assume that the thickness of 
the shell wall is small and that the radius of curvatur 
of its surface is everywhere large compared with the 
wave length. The latter hypothesis allows the shell to 
be approximated as plane in establishing the interaction 
of the surface wave with the external air, in consequence 
of which we may assume a plane wave propagating out- 
ward from the surface (this approximation is examined 
in Appendix A). 

It is convenient to express the equation of motion of 
the shell in a coordinate system (x1, i, 21, 4) moving 
with the flight velocity and that of the surrounding air 
in a system (x, y, z, /) that is at rest relative to the air at 
infinity. In the first system z; = O defines the undis- 
turbed position of the shell, x, and y; are measured in 
the surface of the shell, and the flight velocity Ll’ is 
directed along the negative x-axis; the two systems 
then are connected by the Gallilean transformation (see 


Fig. 1) 
m=x+Ut M=y 4 = 3 A= Cis 
x=x,—- Ut y=a 3= 2% t= (3.1b) 


We emphasize that partial differentiation with respect 
to time is not invariant under this transformation; on 
the other hand, vy, and 2; may be replaced by y and 2 
throughout the subsequent analysis. 

The linearized equation of motion for small displace- 
ments (¢) of the shell relative to zs = 0 may be posed in 
the form 


— Db, (3.2) 


Le + m(07¢/dt?) = p; 


where L is an operator (which may involve differentia- 
tion with respect to the space coordinates but not time) 
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SUPERSONIC FLUTTER 
such that L¢ is the transverse load associated with the 
deformation ¢, m is the mass per unit area, and p; and 

are the internal and external pressures, respectively. 
The question of how small ¢ must be to justify the neg- 
lect of nonlinear terms in the operator L has been ex- 
amined by Reissner.* He concludes that if bending is 
must be small compared with the wall 


umportant |¢ 
thickness, but if pressurization effects dominate those 
of bending |¢) need be small only relative to the static 
displacement of the shell. The neglect of nonlinear 
terms in the following treatment of wave motion in the 
external and internal fluids also requires that ¢ be 
small compared with the wave length and that the 
changes in density and pressure be small. 

A general form for a plane wave in the Cartesian co- 
ordinates x, y, 3 1s given by (the real part of 


' } 


n lL exp ) (277 A) (al v sill # cos a 


( Oo” 


ysin@sina —scos@)y (5.5 


where 7 denotes the displacement normal to the wave- 
front, -1 the complex amplitude, \ the wave length in 
the fixed (x, ¢) coordinate system, @ the (polar) angle be- 
tween 7 and the z-axis, a the (azimuthal) angle between 
the x-axis and the projection of » on s = O (see Fig. 2), 
and a the speed of sound in air. The angle @ must be 


real and smaller than 7 2 


in magnitude insofar as 
represents a wave propagating outward from s = 0; 
however, Eq. (3.5) also remains valid (as a solution to 
the wave equation) for properly chosen complex values 
of 6 [see Eq. (3.9b) ]. 

Che formulation of our problem requires the matching 
of the transverse component of the plane wave defined 
by Eq. (3.3) with the transverse wave ¢ and the calcu- 
lation of the resulting acoustic pressure. This condition 
may be imposed at the plane s = O (rather than s = ¢) 
in virtue of the assumption of small displacements, 


whence 


¢C = (#2 cos 6 -1 cos 6 exp 


‘ . ° ° ' » 
} (277 A)(at — x sin @cosa — ysin@sina); (3.4) 


Alternatively, we may say that the transverse wave of 
Eq. (3.4) would produce the plane wave of Eq. (3.5) in 
>. 
length (A,) of the transverse wave of Eq. (3.4) 


The amplitude (.1,), wave speed (c), and wave 
in the 


two-dimensional coordinate system (x, y, f) are given by 


A, = A cos 0 3.4) 
ri = acsc @ (3.6) 
A. = Acsc#@ = (c/a)aA (3.7) 


and, in terms of these parameters, we have 


. é ° ° i ‘ 
¢ = A,exp } (277 Ai)(ct — x cosa — ysina)y (3.8) 
cos @ = Vl — (a/c)* (3.9a) 


We note that \; is, by hypothesis, a positive real num- 
ber throughout the remainder of this analysis. This 
imposes no essential restrictions on the end results. 

If c¢ 
surface wave, similar to the refracted disturbance follow- 


< a, the disturbance of Eq. (3.5) becomes a 


OF 


ACYLINDRICAL SHELL LOY 


ing total reflection at an interface ;!'’ it then follows that 
in satisfaction of the requirement that ” must bs 


> 0 


bounded throughout the half space s 
cos 6 ulacv i c/a)s ( 3.9b 


3.9b 


radical not only for real c but also as an analytic func 


Eqs. (5.9a) and define the principal value of the 


tion of ¢ in a complex plane cut from « tato +a 


exp (iT 2 


The acoustic pressure asst ciated with the plane wave 
n must be equal to the rate of change of momentum 


across the wave front and is given by (the linearized 


approximation 


Pp. p,a(On Ol 2r1/Xr 3.10 


p,ad-u 


where p, denotes the density of the undisturbed ai 


The pressure on z = O, obtained by expressing m in 


terms of ¢ from Eq. (3.4) and introducing c¢, A;, and cos 4 


from Eqs. (3.6), (3.7), (3.9a), and (3.9b), then is given 
by 

p 2 ripyacd,—'{1 a ¢ ¢ 3.1la 

2 paC-Ay l — (c/a ~* 3.11b 

We now return to the equation of motion, Eq. (5.2), 

in the local coordinate system moving with the shell. 


Ihe surface wave ¢ in this system, obtained by substi 


tuting x and ¢ from Eq. (3.1b) in Eq. (3.8), is given by 


l(c +l 


— 7 


. )/s 
( I, exp } 21 Ay cos ajl 


i ‘ 
COS a ¥ Sil a}, ee 


We find it convenient, in the subsequent analysis, to 
focus attention on the axial phase velocity relative to 


that speed at which an observer 


the shell—that is, 
would have to move along the positive x,-axis in order 


to view the disturbance ¢ as a standing wave Denot 

ing this speed by |’, we have from Eq. (3.12 

¢ = 1, exp ; —(271 Ay ie Vt,) cosa ay sin a { a.13 
( 4 a= l COS Qa 3.14 


The transverse acceleration corresponding to Eq. (3.15 


is then given by 


o°¢ Of” —((2eV Ay |°¢ 3.15 


COS @ 


The pressure /; resulting from the disturbance in the 
internal fluid may be determined by substituting c in 
Eq. (3.11) from Eq. (3.14), replacing ¢ by —¢£ (since —¢ 
is the motion info the internal fluid), and setting / 0) 
(since the internal fluid is at rest relative to the shell in 
its equilibrium position). The result may be further 
simplified if the internal fluid is incompressible 
» ; since the effects of the internal fluid are not likely to 


’ 


1.e.,,a= 


be important unless it is a liquid, we incorporate this 


approximation. We then have from Eq. (3.11b 


2rp;( V cos a)*A\,~'¢ 5.16 


= 

We remark that p; represents the virtual inertia force of 
the internal fluid. 

The structural operator L may be related to the axial 

wave speed 1) that the transverse wave ¢ would have in 
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the absence of the aerodynamic forces. Substituting the speeds l’ + a sec a—equivalent to speeds ¢c = +a 
Eqs. (3.15) and (3.16) in Eq. (5.2) and setting p, = 0, relative to the air at rest—correspond to plane waves of 
we obtain sound having normals parallel to the shell and not in- 
, teracting therewith (we may imagine that u = O implies 
Le = [m + (pi 27) (27 Vo cos a \y)76 0 (3.1% ” ‘ ez i 
m = o, rather than p, = O, for these waves). If we 


We emphasize that |) is an eigenvalue of the partial 
differential equation |Eq. (3.17) | and will exhibit a de 
pendence on A, and a that must be determined by the re 
quirement that Eq. (3.13) shall be a solution to Eq. 
(3.17) when I” = I. Wealso note that |) will be real 
only in the absence of structural hysteresis. 
Substituting Eq. (5.11b) and Eqs. (3.15) through 
(3.17) in Eq. (3.2), we obtain 
[7 + (pjAy 27) |(2m% cos a X,)*( Vo? — V?) = 


” 


2m pac” Ai Vv 1 — (c/a)- (3.1S 


Introducing the mass parameter (which differs from that 


of reference | 


KM = pad: (pit1 + 27m) 3.19 


and substituting c from Eq. (3.14), we may rewrite Eq. 
(3.18) in the form 


V2 + w(V — U)241 — ((V — U) asec al?t 7”? = Ve’ 


It remains to determine |) asa function of A; and a from 
Eq. (3.17) and to solve Eq. (3.20) explicitly for I’, in- 
terpreting the radical as a single-valued function of I’ 
by drawing branch cuts from |’ = Ul’ + a sec a to o 
exp (/m 2), in accordance with Eqs. (3.9a) and (3.9b 
The nature of the panel motion then may be inferred 
from Eq. (3.13); in particular, the motion will be un- 
stable if the imaginary part of I” is negative. 

The generalization of Eq. (3.20) for nonplane waves 
and its solution for /ong wave lengths (compared with 


the radius) is discussed in Appendix B. 


(4) THe STaBILIty EQUATION 


The stability equation |Eq. (3.20) ] may be reduced to 
the corresponding two-dimensional equation of ref- 
erence 1, of which a Cauchy-Nyquist analysis revealed 
that for |) rea/ it had either zero or one root with nega- 
tive imaginary part as U was less or greater than a 
critical speed Ll’... The results of that analysis are 
directly applicable to the present problem if .1/ is re- 
placed by .\f cos a, but we now wish to include the 
effects of structural hysteresis (Il) complex). We also 
shall impose the approximation 4 < | a priori, thereby 
greatly simplifying the analysis without appreciably 
limiting the practical application of the results. 

We begin by noting that Eq. (3.20) may be rewritten 


in the form 
(V — Vo)(V + Vo)V V — (U — asec a) X 


V(U + aseca) - V+ pa seca(V — UU)? = 0 
(4.1) 


If we set u» = Oa priori, this equation has the primitive 
roots V = +V,and V = U + asec a; the speeds + Vp 
correspond to transverse waves on the shell alone while 


make yw approach zero through positive values, however, 
we find that the cuts from 1° = lL’ + asec ato © exp 
(17 2) restrict the possible wave motions as follows 


The transverse wave that travels upstream relative to 


the shell with speed approximating |\(l7 ~ Vo}) is 
damped by radiation if |Ii, > a sec a U’ [see Fig 
3(a), (b), (d), (e), (f)]; if | Vo) < asec a — LU [see Fig 


3(c)| no radiation occurs, and the only effect of the 
aerodynamic forces is a slight decrease in speed (access 
of inertia Structural hysteresis always serves to 
damp this wave. 

The transverse wave that travels downstream relat. 

to the shell with speed approximating I, is damped by 
radiation if |1% > U + asec a@ [see Fig. 3(a), (d)], is 
undamped but slowed down by aerodynamic inertia if 
U —aseca< |Vo| < U + asec a [see Fig. 3(b), (c 
(e)], and acquires energy from the air stream if |’ 
l’ — a sec a@ [see Fig. 3(f)]. This last possibility exists 
only for supersonic speeds. Structural hysteresis also 
tends to damp this wave, but, as we shall see subse 
quently, it may not suffice to eliminate instability if 
approximates Ll’ — a sec a. 

The wave that travels downstream relative to the air 
at rest with speed approximating a acquires the charac- 
ter of a surface wave (cos @ imaginary) in the presence 
of the shell if | 1) > U' + a sec @ [see Fig. 3(a), (d) 
and is damped by structural hysteresis; if | I) < U' + 
a sec @ [see Fig. 3(b), (c), (e), (f) | this wave is no longer 
dynamically possible. 

The wave that travels upstream relative to the air at 
rest with speed approximating a also acquires the charac- 
ter of a surface wave in the presence of the shell if 

Vo| > |U — a sec a| [see Fig. 3(a), (b), (d), (e)]; 
structural hysteresis tends to damp this wave only if 
U <a sec a@ [see Fig. 3(a),(b) |, however, and actually 
leads to a weak instability* if l’ > a sec a@ |see Fig. 3(d), 
(e) |—again a possibility only in supersonic flow. This 
wave also is dynamically impossible if | ly) < |U’ — a 


sec a| [see Fig. 3(c), (f) ]. 

These results are illustrated in Fig. 3(a)—(f) for the 
possible ordered sequences of the primitive roots; the 
positions of the roots for « = 0+ and real I (no struc 
tural hysteresis) are shown as open circles, those for 
small positive u and real |) as dotted circles, and those 
for small positive « and small structural hysteresis 
11,2 = | Vy|? exp [7g sgn (RI I’)]{ by solid circles. 

We emphasize that these conclusions are predicated on 
the assumption of small u and are not uniformly valid 


with respect to the proximity of either + 1) with U’ + a 


* The damping ratio for this instability is given by 


6 = —gyu*(|Vo\ cos a/a)*( MM cos a — 1 (|Vo| cos a/a)? — 
(M cos a — 1)?]~3 (4.1’ 
This is quite small unless Vy approximates Ll’ — a see a, in which 


event the analysis of Eq. (4.2) et seq. is applicable 
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sec a. Nevertheless, they do indicate that the most 
serious instability for supersonic speeds will be asso- 
ciated with the confluence of the two primitive roots | 
and L’ 
hood to a more refined analysis. 

if both | and V are assumed to approximate ( a 
sec a, We may (approximately) factor out (V7 + I) and 
(U + asec a — V)'? in Eq. (4.1) to obtain the ap 


asec a, and we therefore subject this neighbor- 


proximation 


a sec a) | , 


(V lV, )[V (l 
2~*/"0(a sec a)” "(U — asec a)! O (4.2) 


We shall assume, in the remainder of this analysis, that 
1) can achieve the value 1’ — asec a and that yu 1s suf 
ficiently small to justify the neglect of u’/°(.17 cos a — 


1) 1° relative to yu” * lef. Eq. (4.3 


2 
3 


(Jf cosa — 1)” 
below |. 

We now seek to establish more precisely the condi- 
tions under which Eq. (4.2) has a root (or roots) with 
negative imaginary part and to calculate the maximum 
magnitude of its phase angle—t.e., the maximum nega- 
tive damping ratio. The analysis is expedited by the 
scale transformations 


|’ (U’ — asec a)[l + (1/2)u(u Dyers ye 
(M cos a — 1)~**], —(32/2) 
argu <7/2 (4.3a 
Vy = (U — asec a) [1 + (3/2)(E + in) (u 2)?" X 
(M cos a — 1)~° " (4.3b 


where v is a dimensionless complex variable that serves 
to expand the neighborhood of the branch point |7 = 
lL’ — a sec a, & measures the proximity of | 1) to this 
point, and 7 measures the structural hysteresis. If g de- 
notes the usual structural damping ratio [such that 
Vo? = | Vo7/(1 + zg), g< 1], we have 

n = (1/3)g(u/2)~* (Mcos a — 1)”* (4-4) 


5S 


Substituting Eqs. (4.3a) and (4.5b) in Eq. (4.2), we 


obtain 


3/2 ‘ & ° 2 ‘ 
wiv) = v — 3(& + inv! — 2 = 0, 


where w is analytic in the cut v-plane. 

A Cauchy-Nyquist analysis of Eq. (4.5) reveals that 
if 7 = O, w(v) has O or 1 zeros in Im v < Oas&> 1 or 
& < 1, respectively, as shown in Fig. 4(a)—(d)._ We may 
compare this more refined approximation to the sta- 


bility boundary in the absence of structural hysteresis 


viz., £ = |——-with the rough approximation — = 0 im 
plied by the dotted circles in Fig. 3 [especially 3(e) and 
(f) |. 


The revision of the Cauchy-Nyquist diagrams re- 
quired by structural hysteresis is shown in Fig. 4(e) 
(h); this analysis leads to the conclusion that w(v) has 
0 or 1 zeros in Im v < Oas & < —(4/27/n’) or E> — 
(4/27n"), respectively. We infer from this result that 
structural hysteresis is capable of preventing instability 


in the neighborhood of 1’ = J) only if 1 is sufficiently 
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small relative to U — asec a [ef. Fig. 3(f)]. If, on the 
other hand, |) is sufficiently close to l’ — a see a strur 

tural hysteresis cannot prevent instability, while if | 

is sufficiently large compared with lL’ — a sec a (£ > | 

instability actually must be charged to structural hys 
teresis in the sense that the motion would be stable in 
its absence [cef. Fig. 3(d) and (e We remark that this 
rather striking result does not imply that energy can be 
produced by hysteresis per se, but rather that the phase 
shift associated with hysteresis may allow more energy 
per cycle to be extracted from the air stream than is 
dissipated in the structure. 

Having established the conditions under which w(: 
has a root with negative imaginary part, we now seek 
to establish the maximum value of the latter quantity, 
corresponding to maximum growth per cycle of the un- 
stable motion. We first separate the real and imaginary 
parts of v, say Band —y, as functions of the independent 
variables — and y —viz., 


v = B(E, n) — t7(E, (4.6 


Substituting Eq. (4.6) in Eq. (4.5), differentiating par 
tially with respect to £, and equating the real and 
imaginary parts of the result separately to zero yields* 


(Oy /0£), = O (4.7 


wtr 


¥ +76 = 0 if 


This, together with the requirement that both n and 4 
be positive, requires the phases of & + im and 6 — 7y to 
differ by +7; under this constraint Eq. (4.5) admits a 
solution in the cut plane only if the phases of these two 


complex numbers are taken as 75 and —27 3, re- 
spectively, and we write 

(4.Sa 
(4.Sb 
Introducing these expressions in Eq. (4.5) yields 
3/2 1/2 _ 9-1/2, 93/4 _ 9 (4.9 


ws 9 
Y — ony — 2 °.) = 


The only positive root to Eq. (4.9) is found to be 


c = cos (7/6) 


We now return to the original variables of Eq. (4.2 
Imposing the constraint (4.Sa) on Eq. (4.3b) and 
eliminating » through Eq. (4.4) yields 

Vo(A, a)| = (U' — asec a)[1 + O(g) | (4.11 
while substituting v in Eq. (4.3a) from Eqs. (4.8b) and 
(4.10) and evaluating the damping ratio yields 


*We recall that the original independent variables, with re 
spect to which instability should be maximized, are the wave 
length A; and the propagation angle a; within the approximation 
u-/3 (M cos a — 1)~°/3€1, the desired end may be accomplished 
by maximizing y with respect to and using the resulting value of 
— to establish a constraint between @ and \;. The explicit repr« 
sentation of this constraint and the completion of the maximizing 


process require Vo to be exhibited as a function of @ and \ 
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6 Im | RIV (4.12a) 
1 2)y(u 2)7? (cos a Pa 
Olu! 7 (AM cos a \-F 1.12b) 
g 3 2-*/F (VG? + (1 2)(g 3 6°27 PF 
2-2/3 14/ 6% + (1/2)(¢/3)? — &°*7 7°" 1.12¢ 
where 6) is the negative damping ratio when g 0 
VIZ 
i V3 4)u°? (i cos a L13a 
V3 4) [pdr (pr + 2am) |? ?% Xx 
A COS Q l _ 1.13b) 


It remains to maximize —6 with respect to A; and a, 


subject to the constraint (4.11), but this cannot be 
done until the explicit dependence of V7 
established. We remark that, for fixed 69, 6 is a nono 


tonically increasing function of g, so that it suffices to 


on A; and a@ is 


maximize 69, after which 6 can be calculated from Eq. 
+. 1 Ze). 


Che ratio 6,6) is plotted in Fig. 5 as a function of 


6b (3 r[V1i + (1 2)(¢ 36)? + 17? 
0 30 | 2 | 
(V1 + (1 2)(g 36)? — 1 |” 
t. 14 
2 
We note the following approximations: 
6 —59 + (1,3)2 + O(g? 6 b.15 
6 = —4g-%)? + O(g-5°) L16a 
—(3V 3 16)(u g)? (Wl cos a — 1) + 
O[utg— (VW cos a — 1)~"] (4.16b 


It follows from Eq. (4.15) that structural damping in 
sufficiently small amounts is two thirds effective (rela 
tive to the value 6 = g/2 that would result in the com- 
plete absence of aerodynamic forces) ; in large amounts 
Eq. (4.16) indicates the best that can be achieved is a 
weak instability of order g~*. (By hypothesis g always 
must be small relative to unity, but it may be large rela- 
tive to 69; if both g and 69 were sufficiently small, insta- 
bility might be prevented by dissipative effects—vis- 
cous, thermal, or relaxation—in the air.) We also note 
that the condition & < —(4 27m") implies that struc- 
tural damping will stabilize the motion when 


ly < (U — asec a)[1l — 


(1/2)p?g—*(. 11 cosa — 1)*] (4.17) 


The foregoing analysis could be extended to cover the 
neighborhood of the pseudosonic point 7 cos a = 1, 
where the transformation Eq. (4.3) requires modifica- 
tion. The first approximation to the stability bound- 
ary—i.e., the counterpart of 1) = U' — a sec a—is 
found to be! 


2) yu? > 4 O(/M cosa —1}) (4.18) 
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Fic. 5. The factor by which the negative damping ratio do, 
calculated on the assumption of zero structural damping, is 
reduced by finite structural damping; the dashed line gives the 
reduction that would occur if the structural and aerodynamic 
damping factors could be calculated independently and then 
added directly 


- . . 1 2 - 
1.5) is a polynomial (in ¢ of 


The counterpart of Eq. 
fifth order, which would require numerical solution for 
the explicit calculation of the negative damping factor 

Again we emphasize that all of the results deduced 
from Eq. (4.2) et seg. assume that I’) approximates [ 
asec a and that wis small. If uw is small and the mini 


mum value of I’) exceeds a by an appreciable 
amount, only the very weak instability of that primitive 
mode for which |’ approximates l’ — a sec a |see Fig 
3(d) and (e) and footnote Eq. (4.1’)] can occur. If u is 
large, on the other hand, the critical value of Ip is / 
(see reference 1), but then the plane wave approxima 
tion may break down for a nonplane structure (see 


Appendixes A and B). 


(5) WAVE SPEED FOR CYLINDRICAL SHELL 


The types of waves that may be propagated along a 
cylindrical shell may be classified in a general way as 
longitudinal, tangential, or transverse as the dominant 
motion is axial, circumferential, or radial. The last 
class of waves is the one of principal interest in the 
present investigation, both because of the intrinsically 
smaller wave speeds and because only transverse motion 
gives rise to aerodynamic forces. 

The equations governing transverse wave propaga 
tion might be deduced from the general equations of 
vibratory motion established by Rayleigh! and 
Love,'* but we remark that, since frequency varies like 
\~! |’, the wave lengths of minimum wave speed and 
minimum frequency must be expected to be quite dif- 
ferent. For this reason, the anticipated deflection pat- 
terns resemble more closely those calculated (on small 
deflection theory) for axial buckling of a cylindrical 
shell, the axial load being analogous to the inertia 
force.* It is known" that such buckling characteristi- 
differ 


markedly from those predicted on small deflection theory; how 


*We note that observed buckling patterns usually 


ever, We anticipate that nonlinear structural effects will be rather 
less important in the dynamical problem Love's re 


at least in the determination of sta 


ace, €.2.. 
marks, reference 12, p. 553 
bility criteria. 
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cally occurs with a wave length that is sufficiently small 
to make those extensional and flexural stresses asso- 
ciated with transverse deflections of approximately 
equal importance but sufficiently large to render ex- 
tensional stresses dominant in the axial and circum- 
ferential directions. Motivated by this result, we shall 
include only direct extensional stresses in the axial and 
circumferential equations of motion, noting that the 
accuracy of the implied approximations may be estab- 
lished a posteriori. We also shall neglect the axial and 
circumferential accelerations (u,, and v,,) and the rota- 
tory inertia of the shell wall. 

The equations of equilibrium of a_ thin-walled, 
evlindrical shell are, on the basis of the foregoing ap 


proximations, |! 


rT) + 1 — pv) 2)u,, 4 1+ pv) 2? 
(vy R)w, = 0 d.la 
(1 + v)/2)uy + (1 — v)/2]e T Uy ~ 
1 R)w 0 1b 
(kh R) (A - 1 —Yur — Vy T 
R-'w) — [Nywr: + Nx - 
Diw + 2w ryy TT Wyy a = Pp 5.le 


where “, 7, and w denote the axial, circumferential, and 
normal (radially inward) displacements of an element 
of the wall, « and y the axial and circumferential co- 
ordinates (v need not be distinguished from .; in deter 
mining |’), & the radius of the cross section, / the wall 
thickness, v Poisson's ratio, & Young's modulus, .V, and 
NV, the axial and circumferential preloads per unit 
length, p the normal load per unit area, and D the bend- 
ing stiffness. The terms in parentheses, brackets, and 
braces on the left-hand side of Eq. (5.1e 
effects of extension, preload, and bending, respectively. 


We seek a reduction of Eqs. (5.1a, b,c) to [ef. Eq. (3.2 


Lw = P ie 


represent the 


for deflections having the spacewise dependence cf. 


Eq. (3.13 


9 ay » , : ( ae 
U,v,;wW~ exp} —(2at/Ay) (xX COS a + Y SIN a) |} a) 
Substituting the displacements of Eq. (5.3) in Eqs. 

5.la, b, ¢) yields algebraic equations in u, v7, and w, 


solving these for w and comparing the result to Eq. 
(5.2), we obtain 
b= Eh R?*) cos! a + (CN, cos? a + N, sin? a) X 


25 Ay 5 D(27r Ay : (5.4 


The corresponding wave speed, as defined by Eq. (3.17 
then is given by 

V (Ay, a = 

Ith Ay cosa 2rR)?+ (N, 7 N, tan°a) + D(2r A; COS @) 


m + (pj, 27 


If we substitute 


a;? = Eh/m = E/p, (5.6 
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D = Eh 12(1 — pv 5.8 
in Eq. (5.5), where a; is the velocity of sound im the 
shell material, p; its density, and e, and e, measures ot} 


the initial strains (€, and «, would be the true strains 


only if y = O), we obtain the alternative expression 
| (A , a = ad," 1] + p \y 2arph x 
) (A; cos a 2rk)? + (e, + €, tan? a) 4 
11/1201 — pv 2arh rd; cos a)*; 9.9 


which is especially appropriate when the shell inertia 
dominates that of the internal fluid (p,A, 2rph 
The result of Eq. (5.5) may be generalized for appli 
cation to an aeolotropic shell by introducing nonequal 
effective values of the wall thickness in axial extension 
h,), circumferential extension (/,), and shear (/,,) and 
of flexural and torsional rigidities for an axial element 


/, and $s 


The transverse 


(J, and /,,) and a circumferential element 
all in the notation of Timoshenko.! 
wave speed, assuming the wave length to be /arge com 
pared with the spacing of stiffening elements (so that 
their contributions to /,,..., /,, may be averaged over 
a wave length) but small compared with the radius and 


length of the shell, is given by? 


i Ai, a) = My, tT \p \y 27 ,E cost a h 
2 sin? a cos? a h,,) + (sintah x 
(\, cos a 2rR)? + CV, + NV, tan? a) 4 

E|I, costa + U + | sin” @ Cos” a 


f sin? a 2r dX; COS a } ().10 


We emphasize that the results just obtained for |) are 
that are smal] compared 
of the shell: 


valid only for wave lengths (A, 
with both the radius (/) and length 
the neglect of axial and circumferential bending and 
inertial forces is valid only if A R, while the neglect ot 
the end constraints requires \; < It is perhaps well 
to recall, on the other hand, that we have implicitly as 
are valid only for 


We shall find 


that these inequalities usually will be satisfied in the 


sumed \,; > h, since the Eqs. (5.1 


thin shells and neglect rotatory inertia. 
neighborhood of minimum (with respect to A; and a 
wave speed, but they must be borne in mind in interpret 
ing the limiting trends of the results for very large o 


very small wave lengths. 


(6) Empty ISOTROPIC SHELI 
We consider in this section the more specific con 
figuration of an isotropic shell that is not prestressed 
We also shall neglect the inertia of the internal fluid, on 


the assumption that 
pitt < 2arph 6.1 


as usually would be true for a gas (but not a liquid 


The wave speed given by Eq. (5.9) then reads 


aes id ? ' 
Vu?(A1, @) = a7} (A, cos a@ 27R - 


(1 12(1 — v?)](2rh \, cos a)?} (6.2) 
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SUPERSONIC FLUTTER 
« + a > 
i , P 
\ a ee 7 
\ 
\ 
\ 
‘ 
= Y 
‘ ral 
i. 
FI 6 The variation of axial wave speed for an empty 
jsotropic evlindrical shell; @ is the angle between the wave nor 


n) and the x-axis, and A, is the wave length in the direction 


f propagation Set Eqs 6.3) to (6.5 


The wave speed of Eq. (6.2) depends on the single 
variable \; cos @ (rather than A; and a independently) in 
a very symmetrical way that may be emphasized by re 


writing it in the form 


Vo" 1 2)V*?[(A; cos @ A,*)? + (A,* A; cos a)?} (6.3 


where |°*, the mini 1 value of I), 1s given by 
! ]*, the minimun lue of | \ by 


; 2 1/47) 1/2 
V* = [301 — pv h R)° “a 6.4 


and \,*, the value of A; cos @ at which this minimum 


occurs, by 


i wa 14 12 a 
Prl12(1 — v*)] “(Rh 6.5 


\i* 


rhis result is plotted in Fig. 6. We observe that ex- 
tensional and flexural deformation make equal contribu- 
tions to the wave speed at the minimum point; more 
over, to any prescribed value of the wave speed greater 
than |’* there correspond two values of A; cos a having 
the geometric mean \,*, with extension or flexure domi 
nating as A; cos a exceeds or is inferior to A,;*. This be 
havior of the wave speed is analogous to that found for 
surface waves on water,'® where gravity and surface 
tension appear as the counterparts of extension and 
flexure, respectively; however, in that problem the ex- 
ponents of (A A,* 1, rather than +2 as in 
Eq. (6.3). We also remark that Eq 
the wave length in the neighborhood of the minimum 


would be + 
6.5) implies that 


wave speed is of the order of the geometric mean of the 
radius and thickness of the shell, whence the hypotheses 
\, < Rand \,; > fh both will be satisfied for thin shells 
h<« R) 

rhe flutter speed—i.e., the minimum flight speed re- 
is given in first approxi- 
mation by the minimum value of 1(A,, a) + 
Regarding A; cos a and a as independent variables, we 
A” 


quired to produce instability 


a sec a. 


find that this minimum occurs at a = 0 and \; = 
and has the value 


U* = V* +a 6.6 


Eliminating I’* between Eqs. (6.4) and (6.6), we find 
that the shell wall thickness radius ratio required to 
prevent flutter at a flight speed U is given by 


h* R) = V30 — v*)[(U —a) a)’, U>a (6.7) 


This approximation is too crude near Ll’ = a, but it 
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suffices if faired into the end point (at / a) deduced 
from Eq. (4.18), which, in conjunction with Eqs. (6.4 
and (6.5), vields 


h* R) = 1.29(1 v2)" “(aa Pa P 

at l a os 
The numerical values of (h* FR 
(6.8) are plotted vs. Mach Number in Fig. 7 for duralu 
an altitude of 35,000 ft 


given by Eqs. (6.7) and 


min and steel at sea level or 
We remark that Eq. (6.7 
velocities a and d,, is relatively insensitive to changes in 
1,000,000 


, depending only on the sonic 


altitude, especially in the troposphere 35.000 
ft.), where a is approximately constant 


If the ratio # Kis inferior to the critical value (/* Kk 

flutter is predicted by the analysis of Section 4, and we 
may determine the critical negative damping ratio by 
maximizing 6), as given by Eq. (4.15 
\; and a, the constraint 


effect of structural hysteresis on the damping ratio then 


with respect to 


subject to +1] The 


may be evaluated from Eq. (4.14) or Fig. 5 
If a factor (cos a)” ” is absorbed in the parameter u 
in Eq. (4.13) and the remaining dependence on cos a 


removed through the constraint [rewriting Eq. (4.11 


COS a a l | 6.9 
we obtain 
6 V3 4)(p,A; cos a@ 2rph)~ Vo /a x 
l - Vi a} 2 l > | T a 6.10 


which, referring to Eq. (6.2), is a function of the single 
We find that a satisfactory approxi 
to the 


variable \; cos a 


mation (less than 7 per cent error maximum 


value of Eq. (6.10) as a function of A; cos a may be ob 
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Fic. 8. The critical damping ratio 69 tor an empty isotropic 
cylindrical shell; see Fig. 5 for the effect of structural damp 
ing 


tained by assuming it to coincide with the minimum of 


lo; substituting \; cos a = \\* and I) = I* from Eqs. 
(6.5) and (6.4) then yields 


&5* = (3/4)°/4(1 — v?)'/*4(2a,/a)*/*(p1/p,)'”? X 
[(pal?R/2Eh)? — (12)7”4(1 — v2) 4 X 
(pa/p)?)"® U>V*+a (6.11) 


The corresponding numerical values are plotted in Fig. S 
for duralumin and steel at three altitudes. 

The frequency at which flutter occurs is, by reference 
to Eq. (3.15), 


f=r"Vecosea (6.12) 


The frequency corresponding to Eq. (6.11) is obtained 
by letting \; = \1* sec aand V = J* in Eq. (6.12) and 
substituting cos a from Eq. (6.9), with the result 


~- V*)F> 
(4/2a, 24R)M 


f* = (V 2a,/27R) [a (U 
(6.13) 
The frequency at the flutter speed (LU = I’* + a) is, 
therefore, inversely proportional to the time that would 
be taken by sound to travel a distance equal to the cir- 
cumference of the cylinder /m the material itself and falls 
off approximately as the inverse square of the Mach 
Number at higiter speeds. We find, for example, /* = 
750 eps for a steel or duralumin fuselage 10 ft in diameter 
at the flutter speed——viz., U = V* +a. 


APPENDIX A—-AERODYNAMIC FORCES 


We investigate here the accuracy of the plane wave 
approximation for a circular cylinder. The transverse 
wave motion along the surface of the cylinder given by 


Eq. (3.8) is, in cylindrical polar coordinates (x, 7, ¢), 
(A-1) 


¢ = A,exp } (2x7 Mi) [cet — x cos a — Resin a} } 


We require a solution to the wave equation 


AERONAUTICAL 
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av-*d = on (A-2 
subject to the boundary condition 


Oo = 6: = (27Ic Ay ¢ r R A-3) 


- 


together with the radiation and finiteness conditions a 
infinity. 

We first remark that if ¢ is to be continuous in ¢, a 
and A; must be constrained according to 


2rR sin a’; = n (A-4 


where 7 is an integer; however, insofar as \) < RK we may 
assume a and ), to be distributed over a continuous 
A solution to Eq. (A-1) 


(A-35) and the conditions at infinity is given by 


spectrum. that satisfies Eq. 


@ = (2m1Rec/d,) [K,(2r/R)/2K,'(2) |f (A-5 
2rR : 
s= V cos- a — (c/a (A-6 
\y 


provided that z is required to have a positive real part 
(finiteness condition) and a non-negative imaginary part 
(radiation condition); A, is a modified Bessel function of 

the second kind. 
The perturbation pressure on the cylinder is given by 
Pa = —p,9,; (A..7) 


Substituting Eq. (A-5) in Eq. (A-7) and introducing the 
function 


y,(s) = —K,(2) /2K,"(z (A-Sa) 
= 2?-1(K.44(s) + K,-a(z 1K ,(z) (A-S8b) 

we obtain 
Pa = —Pal27rRe’rj)*¥n(s)-(¢/R (A-9 


We may show that the y, functions have the asymp- 
totic approximation (vide infra) 


vni(z) ~ (mn? + 32)7 1? (A-10) 

while for small z we find 
(2) = In (2/2) — 0.577 + O(27) (A-11a) 
¥n(2) = (1l/n) + O(2? n=>1 = (A-11b) 


It follows that the approximation (A-10) is exact as 3 
tends to either 0 or © unless x = O; if m = O there 
should be a logarithmic singularity at the origin rather 
than the z~! singularity exhibited by Eq. (A-10). The 
corresponding approximation to the pressure, obtained 
by substituting Eq. (A-10) in Eq. (A-9) and evaluating 
n from Eq. (A-4), reads 


Pa = —(27p,C? Ay) [1 = fC a)*] , ¢ (A-12 


in agreement with Eq. (3.11b). 

It remains to investigate the approximation (A-10) 
at the flutter point. We note first that the function y, 
defined by Eq. (A-Sa) satisfies the Riccati equation 


2(dy,/dz) = (n? + 37)~,? — 1 (A-13 


from which it may be inferred that the validity of Eq 
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A-10) demands 
n? + 32°"? > n+ 1 (A-14 
We have, from Eq. (4.8a), 
 — U) cosa = —a[l — (1 /2)u(u/2)*" X 
(lf cosa — 1)~*%] (A-15 


which, in conjunction with Eqs. (A-4) and (A-6), yields 


3/2 9/3 


n> + 3° = |v 


*(2rR Ay 2) pk [ pil 7 
(A-16 


. ) 
(pA, 27) |]; CW cos a — 1) 


We find, for example, that for the unpressurized shell of 
Section (6) the right-hand side of Eq. (A-16) is O[(p,a= 
pi:)(R/h)”’? sec a], whereas n? is 0[(R/h) sin? a], where 
a is determined by the constraint Eq. (4.11). We con- 
clude that the inequality Eq. (A-14) probably would be 


satisfied adequately at the minimum flutter speed (a = 





()) but might be violated for large a. 


APPENDIX B—LONG WavE-LENGTH INSTABILITY 


Leonard and Hedgepeth” recently have given a panel 
flutter analysis of a circular cylinder in which the plane 
wave approximation is not made; indeed, they imply 
that the most serious instabilities are to be anticipated 
at axial wave lengths—i.e., \; see a—comparable to or 
greater than the radius of the cylinder. This conclusion 
is, in OUr Opinion, open to the following major objec- 
tions: 

|) The waves corresponding to n = 0 (pure dilation 
of cross sections) and m = | (cross sections remain 
circular but execute translation) are rather arbitrarily 
excluded by Leonard and Hedgepeth, ‘‘because neither 
of these two motions involves panel action’’; neverthe- 
less, as we have shown in the preceding analysis, axially 
symmetric motion generally yields the minimum flutter 
speed for an isotropic shell and this at a wave length 
small compared with the radius. 
2) The comparative importance of the long and short 
wave-length stability boundaries for fixed » > 2 (re- 
sulting from the two minima of |’) as a function of the 
axial wave length \; sec a) depends, aside from the 
question of structural constraints such as rings, on the 
relative magnitudes of the respective negative damping 
factors. We shall calculate the long wave-length damp- 
ing factor and show that it is likely to be much smaller 
in magnitude than the short wave-length factor of Eq. 

$.15). 

If the plane wave approximation to the aerodynamic 
pressure is not made, and p; = 0, we may deduce from 
Eq. (A-9) that the required generalization of Eq. (3.20 


1S 


)*y,(Z) == 1” (B-1) 


lf’U>l 
in accordance with Eqs. (3.9a) and (3.9b), write z in the 


1? + (p,R pih)(V — l 


where 3 is given by Eq. (A-6). + a we may, 


form 


Z pV [(U — V) al? — le os 


2rReosa dr, (B-2) 


OF 
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where 7? is a dimensionless, axial wave number in the 
notation of Leonard and Hedgepeth. We require an 
expression for the imaginary part of y, on the assump 
tions that z <1 and» > 0 (since, form = 0, Vy has 
only a single minimum, and the choice between long 
and short wave length need not be made). 

We first remark that y, 


then, if (by hypothesis) the phase angle of 1 is small the 


depends only on 3° and In z 
phase of s (for lL’ > I’ + a) will be approximately — 7 2 
(cf. Eq. (B-2) |, and the imaginary part of y, must arise 
through the logarithmic terms. Introducing the series 


expansion’ for A,, in Eq. (A-Sa), we obtain 


In sz + 


‘2 Ins (B-3 


O(s*, & 


for the leading terms of interest. Substituting z from 
Eq. (B-2), the leading real and imaginary terms are 


given by 
(9 /2)°"} [(l l) a — || 


(B-4 


We proceed further on the assumption that p,R < 
np,h, so that I’) serves as a first approximation to J" [ii 
this approximation were not justified the first approxi 
mation to |’ could be determined by substituting y,, 
Eq. (B-1) then 


n~' in Eq. (B-1) and solving for V’ 
vields 
V = Vojl — (1/2 ((U — Vo) Vo) y¥ 


+ O(p,R ph at 


(p,R ph 
(B-5 


as the next approximation to |’. The approximation to 


the damping ratio inferred from Eqs. (B-4) and (B-5) is 


6 = —[x0"" 2°"(n!)2](p,.R pik) X 
(U — Vo) Vol ll I) a - ]{’ 
a>0O, *<1, U>Ve+ea (B66 


It is evident that the long wave length, negative 
damping factor given by Eq. (B-6) generally would be 
much smaller than the short wave-length factor of Eq 
(4.13). 
tropic shell of duralumin for which h RK = 


We consider as a numerical example an iso 
0.0025, AZ = 

2, and the properties of the air correspond to sea level. 
Fig. 7 yields a minimum Mach Number of 1.5 for in 

stability, while Fig. 8 yields a negative damping factor 
of 0.08 at @= 2. The result of Eq. (B-6) cannot be 
maximized with respect to the parameters ? and m with- 
out violating the approximations on which it is based 
(which, in itself, implies that maximum instability is 
not to be sought at long wave lengths), but if we set |" 
equal to its minimum value of 0.5a for the specified data 
and use the values nm = 2 and 7 = 0.11 [the smaller root 
of Eq. (17), reference 17] deduced by Leonard and 
Hedgepeth, we obtain —f = 1.9 X 10 It seems un 
likely that an instability as mild as that indicated by 
this last figure could occur in a real fluid. 
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Tables and Curves 


tor 


Detormations and 


Stresses in Circular Cylindrical Shells Under 
Localized Loadings' 


JOSEPH KEMPNER,* JAMES SHENG,** ann FREDERICK V. POHLE*** 


Polytechnic Institute of Brooklyn 


SUMMARY 


The results of calculations of displacements, rotations, mem 
rane stresses, and bending stresses corresponding to line load 
ings along generator segments of simply supported thin-walled 
ircular cvlindrical shells are presented in tables and curves 
[The computations are based upon a previously developed solu 


f the Donnell equations and are presented for radius-thick 


ness ratios of 10 to 1,500 and radius-length ratios of 1/5 to 10 
The da 


stresses 1 


determination of the displacements and 
loads 
distributed 


the 


i permit 


the neighborhood of the applied loads. These 


ilong generator consist of a uniformly 


segments 
with respect to a generator midpoint; a 
load, 


symmetric 


radial load, svmmetric 


uniformly distributed radial antisymmetric 


the 


piece W1S¢ 


with respect to same point; and a uniformly 


distributed moment loading. Thus, the loads are statically 


equiy ilent to a concentrated radial force and to concentrated 


ixial and circumferential bending moments applied at the mid 
point of a generator 


Procedures are also described which facilitate the application of 


present results to eylinders under arbitrary radial surface 
ds 
SYMBOLS 
D Eh 12(1 v= 
i = elastic modulus 
= length of evlinder 
M = axial and circumferential internal moment re 
sultants per unit length, respectively 
- internal twisting moment resultants per unit 
length 
V.* = applied axial and circumferential moment r¢ 
sultants, respectively 
, = applied inward radial load resultant 
Ox, Y = arial and circumferential internal transverse 


shear resultants per unit length, respectively 
= mean radius of cvlinder 
wall thickness 
ratios of the axial, circumferential, and inward 
radial displacements to the radius a 
= coordinate distance measured in axial direction 


divided by the radius a (origin at midpoint 


of loaded generator 
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26 = loaded length of generator symmetrically lo 
cated about a ) 

b = Poisson's ratio (taken as 0.3 in calculations 

Ox, Ty, Txe axial, circumferential ind shear membrane 
stresses, respectively 

a a = maximum (extreme-fiber ixial and = circum 
ferential bending stresses, respectively, posi 
tive at inner surface of shell 

¢ = circumferential angular coordinate (origin i 
loaded generator 

= comma under a symbol and followed by either 

x or g denotes partial differentiation with 


respect to the indicated coordinat 


INTRODUCTION 


. CYLINDRICAL SHELLS are used extensively 
in all major industries such as the chemical 
petroleum, electrical, and, more recently, nuclear 


reactor industries. They are also widely applied as 
structural components of submarines, airplanes, guided 
missiles, and rockets. In most applications, loads of 
a localized nature due to gravity and inertia effects or 
to thermal expansion of a connected piping system are 
These loads generally 
The attend 


the 


transmitted to the shell walls. 
consist of force and moment components 


ant local stresses which must be combined with 


basic stresses—e.g., due to internal pressure—can of 


themselves be significant and in many cases require 
careful consideration for the economical design of 


attachments and reinforcements 

Recently, an investigation of the deformations and 
stresses caused by localized radial forces and axial 
and circumferential bending moments applied to thin 
walled circular cylinders of constant wall thickness 
was carried out for the Knolls Atomic Power Labora- 
tory.*-" In this work, the loads were considered to 
be distributed along a segment of a generator, and the 
cylinder dimensions (radius-thickness and_ radius 
length ratios) were assumed to be in a range such that 
the simple Donnell 


brief history of solutions of related problems, including 


equations were applicable. A 
a representative bibliography, can be found in refer 
Recent contributions to the solutions 
P. Bij 


loads were 


ences 2 and 8. 
of such problems have also been made by P. 
laard.''! '* In this latter work, the 


sidered to be applied over a small rectangular area of the 


con 


shell surface, and the equations used in the analysis 
differed somewhat from those of Donnell. 

The purpose of the present paper is to present the 
results of calculations of displacement and stress quan 
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tities caused by the action of localized radial forces and 


L/2 " axial and circumferential moments applied as dis 
L/2 T. tributed line loads along a segment of a generator oj 

we ad . . . . 
5 ages a simply supported thin-walled cylinder. It is show; 


that the data presented can be utilized to approximat 
surface loadings of an arbitrary nature. Except for 
minor refinements indicated in reference 13, th 
a method of calculation is that given in references 2, 3, 4 ' 
and 8. The results are presented in tables and curves 
in a form readily applied by the designer for the est; 
mation of deflections, rotations, membrane stresses 
and bending stresses at and in the neighborhood of th 
(a load resultants. The geometric cylinder parameters 
) radius-thickness and radius-length ratios, respectively 


cover the ranges 10 to 1,500 and 1 ‘5 to 10. ; 


METHOD OF ANALYSIS 


Solutions of the Donnell equations!? were obtained it 
references 2 and S for distributed line loads whos 
resultants are a radial force P*, an axial moment \/,* 
































and a circumferential moment .1/,* (see Figs. 1 and 2 
The loads were applied along a length 26 of the ¢ = 0) [ 
generator as indicated in Fig. 2, where ¢ is the circum 
ferential angular coordinate. They consisted of a uni 
(b) formly distributed radial load, symmetric with respect 
to the generator midpoint (x = 0) [see Fig. 2(a 
;. 1. Sign conventions for coordinates, displacements, m : : Ts eee - 

Fr ign conventions for coordinates, displacements, mem piecewise uniformly distributed radial load, anti 
brane stresses, and internal moment and transverse shear result ; ; : . 
ants symmetric with respect to the x = O plane |see Fig 

2(b)]; and a uniformly distributed moment loading 
[shown vectorially in Fig. 2(c)], symmetric with 
respect tox = 0. 
« 
Pp * cy ee Ce eee ere pe 
lL. me i“ a L . Zoe ; rhe radial force, axial moment, an 1 circumferential 
Ms 2 2 2 moment loadings, respectively, therefore correspond t 
ataais os a distributed symmetric radial load (P* case), de : 
1 (p* aes fined as } 
| Wy"? a” /6) 
ry, P*/25) when x| < (/a) / i 
| X f x i 
. es 2 : a } ° 
(-m%, /é ) 0 when (L2a) > |x| > (6/a)\ | 
a distributed antisymmetric radial load (/,* case 
given as 
(M,* /6?) when 0O<x < @/a) ( 
' 
(a) (b) — (\/,* 6?) when —(6a)< x <0 2) | 
0 when (2, /2a)- > \x| > a)) | 
and a distributed circumferential moment loading 
(.\7,* case) 
be L a rere 
2 Ms 2 "| (.M,* 28) when x) < (6/a) } 
ental P ” 
0 when (L/2a) > |x| > {6/a)\ 


in which L and a, respectively, are the length and mean 
radius of the cylinder, and x is the axial coordinate 





} 
fo 
=, 
$ 
3 





divided by a. 

Formulas for stress and displacement quantities | 

corresponding to the three loading cases are given i 

references 2, 8, and 13 for conditions of simple supports 

at the ends of the cylinder. These conditions imply 

(c) that the radial and circumferential displacements as 
well as the axial membrane stresses and axial bending 











Fic. 2. Applied force and moment loadings. 
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TABLES AND Ct 


moment resultant vanish at x = +(L 2a Con 
venient methods for performing corresponding compu- 
are presented in references 1, 3, 4+, and 13, and 
the results of 


Donnell’s equations are applicable are 


tation 
discussions of conditions under which 
solutions of 


given in references 6, 9, 10, 14, and 15. 


RESULTS AND APPLICATIONS 


) 


[he formulas of references 2, 3, 8, and 13 were used 


to compute deflections w, rotations w,, and w,,, normal 


and normal membrane 


bending stresses o,, and 
and o, for a wide range of geometric shell 
For each of the three loading 
cases Characterized by Eqs. (1) to (3) (see also Fig. 2) 
the ratio of the loaded length to shell length (26 L 
was fixed as (1 20). Because of the large number of 
parameters involved in the analysis of the cylinders, it 
was deemed expedient to so limit the loaded length 
The small 


Ooh, 
stresses ¢ 


parameters (see Fig. 1). 
~) 


rather than any of the shell parameters. 
finite loaded length was chosen so as to bring out the 
salient effects of localized loadings on cylinders in a 
realistic manner; the use of truly concentrated load- 
Q) always introduces infinite bending stresses 
Moreover, it will be shown 


ings (6 
under the applied loads.> 
later that the present results can be used to analyze 
the effect of loadings other than those specifically 
treated. 

Calculations were performed for values of the radius- 
length ratio (a L) of 1.5, 1, 5, and 10. These values 
correspond to moderate to very short length cylinders, 
where the length of a cylinder in practical applications 
is considered to be the distance between adjacent rein- 
or bulkheads. Computations  corre- 
5) were performed for radius- 
100, and 200; and 
of 10, 20, 
Thus, the 


forcing rings 
sponding to (a L) = 1 
thickness ratios (a h) of 10, 20, 50, 
for (a L) = 1, 5, and 10 for values of (a h 
0, 100, 200, 500, SOO, 1,200, and 1,500. 
range of (a h) amply covers that usually encountered 
in shell applications with the exception of very thick- 
walled vessels. It may be noted that when (a L) = 
15, no calculations were performed for (a h) > 200. 
rhe calculations were so limited, since for shells with 
aL) < 1-4, say, it is unlikely that very large (ah) 
values would be permitted unless lateral supports, such 
as rings or bulkheads, were introduced. The resulting 
relatively short effective lengths permit the use of 


larger (a L) values in the cylinder analysis. 


Displacements and Stresses at x = ¢ = 0 


The results of the calculations are presented in Figs. 3 
to 6 and in Tables 2 to 21. Figs. 3 to 6 correspond to 
stress and displacement quantities computed at the 
location of the load resultants, that is, at x = ¢ = 0, 
see Fig. 2). The nondimensional form for each quan 
tity permits the corresponding load case to be readily 
identified. Thus, whenever P* appears in the de- 
nominator, the data given correspond to the radial 


loading (P*) case, etc. 


RVES FOR DE 


KE 


ORMATIONS AND ST 


The displacement quantities plotted may be con 


sidered as flexibility (reciprocal stiffness) quantities. 
Each is essentially the displacement or rotation under 
a load resultant in a direction corresponding to the 
resultant. These quantities should be quite useful 
when dealing with piping or similar problems in which 
the shell affects the amount of load 
Phe 
data show that, in the range of parameters considered, 
the flexibility in the circumferential direction is always 
rhis result, 


the stiffness of 
introduced into the shell and external structure 


greater than that in the axial direction. 
of course, is influenced by the difference in the manner 
of loading in the .1/,* case and .1/,* case 

The remaining data in Figs. 3 to 6 correspond to the 
and membrane 


circumferential bending 


0 for the P* case 


axial and 


stresses at x = ¢ = Because the 
\/J,* loading is introduced as an antisymmetric radial 
load [see Eq. (2) and Fig. 2(b) |, the deflections as well 
as normal bending and membrane stresses for this case 
are identically zero at x = 0; hence, these quantities 
are not indicated in Figs. 3 to 6. Similarly, since the 


\J,.* loading is considered as a distributed moment 


along the ¢ = O generator [see Eq. (3) and Fig. 2(c 
the 
identically zero at ¢ = O, 


stresses can be obtained readily with the aid of Eq. (3 


deflections and normal membrane stresses aré 


and the normal bending 


as follows. At ¢ = O, the internal moment resultant 


JJ, is one-half of the applied distributed moment 


Hence, from Eq. (3 


(V,).-0 = V,* 46 when 4 6 a | 

| 
U.).-0 = 0 when (L 2a > 1X > (6 a)\ 
Since maximum (extreme-fiber) bending stresses are 


(6 h*) times corresponding moment resultants 
result 


simply 
see, for example, p. 40 of reference 16 Eqs. | 


in the following formulas for the inner surface circum 


ferential bending stress at ¢ = 0 


IT ob (M,* a“) \lo=0 = (3/2 Lé anu a/L / 


when vi < (6a . 
») 
|Our (M,* a’)|,=90 = O \ 
when (L/2a) > |x| > (6a 
Furthermore, within the accuracy of the Donnell equa 
tions, 
M, = —(D a w, ~ ww...) } 
¢ ( 
M, = —(D/a) (won + vwyzr) § 
Since w is identically zero at ¢ = O, and hence, (% 0 


9, and hence 


=e yield (M/,),-9 = v(.V, 
the inner surface axial bending stress at ¢ 


Eqs. (6 
0 1s given 


by 


loz (M,* a’ leo=zQ = 


Because the bending stresses given by Eqs. (5) and 
(7) are related in a simple manner to the various param 
eters involved, they are not plotted in Figs. 3 to 6. 
However, it may be noted that, since these stresses are 
simple power functions of (ah), they would plot as 
that the 


straight lines in these Figures. It is seen 
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TABLE 1 
Values of Coordinates Used in Calculations 
8/: TABLE 5 
= Variation of o,/(P*/a?) with (a/h) and x for ¢ = 0; 
1 P) > = 9, x 2 1/40) 
é b/a b/a =(2/cXL/s aie 
——_—_ -——- ———_-— — - _ — % x x 
- — - - e/L=1 
TABLE 2 83 7x10 
Variation of w/(P*/Ea?) with (a/h) and x for ¢ = 0; (6/L = : pe 4 = 
eo . ne 
0 x, X> X3 x, si 
i Te ” a ee = earn a. aa oe weet AO. eames: 
a/h a/L = 1/5 ie ren nee 2 yet . 
3x1 c x 3x10° 
-_ 7 pers eben Sak. Mid S| 7 a3 c 7 7 . 
10} 1.39x102 1.37x10 = -1.32x10% = 1.6x10% = 487x107 1 Ld 1.45% 1.32x1 
: 3 ; me , ae ik 5 
20] 6.55x107 6.4410 6110" 52x10? 2. 29x10* Ria gs . sit oe" 
50 | 5.06x10? 4.95x10? 4.65x103 698x102 —-1.79x10? eaten Mac = tox : 
100 | 2.36x104 2.32n10 = 217x210" = 1. 86x04 = 8. 50x10" 
200 | 1.1110 1.08x10?  1.02x107  8.80x10  9.74x104 
a/L = 1 
20 2494x107 2.93x10* 2.89x10° 2.76x10" 1.66x10° TABLE 6 
50 2. 33x10? 2.31x107 227x102 2.12x10 1.56x10~ Variation of —o./(P*/a?) with (a/h) and x for g¢ = 0; (6 = 
100 l.llx10* 1.10x104 1.07x104 = 9.75 x10? 3285x102 1/40 
500 | 4209x107 4.0210  3.82x107 = 3.30xl0”_—:1.43x109 > 
1500 | 4.74x10° = 4.64x10 4.35x10° —-3..72x10 1.69x10° * 0 xy X> X3 x, 
eels cgi nmatisinninnaliaaipeeneite 
TABLE 3 1 1.34xl0~ 75x10~ 
Variation of ¢,,/(P*/a?) with (a/h) and x for ¢ = 0; (6/L = 2 70x10" >. 116x107 
1/40) 50 3.40x10~ 78x10~ 
= 100] 1.05x10? =: 9.04x10* 
x 0 xy X X x, 200 xl0F 257x107 
a/h = a/L = 1/5 _ — —_— - . 
, 20 865x107 3.59x10~ 
10} 9.64x10* a.95x10t_~—5,a5x10) 125x102) 2.46 a «aut sane 
20} 3.10x10* 2.87107 = 1.47xl0* = 220x104 (9.35 a tat ssl 
50] 1.3710? 1.2610? 7.25x102 ~—s 80x02 ~—S 5. 10x02 soo] 4.190 3.8510" 
100 | 4.07x10? 3.92x103_ = 2.21x103 = 54xl0* = 1. 74x10? 1500) 2.89x102 «2. 51x10 
200 | 1.21x10 1.21x104 7.08x103 = 2413x102 ~S 311x102 -___— 
a/L = 1 
20 | 676x102 6.45x10" ss 5.07x10% Ss -2..75xl0* = -1..47x10° 
50 | 3.56x10? 3.37x107 = 252x107 = 112x103 ~— 138x103 
100 | 1.234104 1.15x10* = 8.17x102 =~ 2. 86x10 1.38x10* TABLE 7 
500 194x010? 1.79x10° 1.09x10° 1.41210 5 .84x102 Variation of w/(P*/Ea?) with (a/h) and ¢ for x = 0; (6/L = 
1500 | 1.14x10 1.08x10°  6.05x105 a.35x10* = 4. 31x04 1/40 
C ) ) % % %, 
TABLE 4 a/h aft, « 1/5 
Variation of og,/(P*/a?) with (a/h) and x for ¢ = 0; (6/L = ; ; —_— 
1/40 10] 1.39x105 1.35x10* 1.25x10° 9 0x10" 4.6910" 
oo eee ee 2 6.55x10 6.27x10 § 63x10 3.93x10* -2.20x10° 
\ x ea - ; 50] 5.06x107 470x103 «-3,93x10?_—S2,.10x10? = -1.11x10 
‘Sree a TT: SI: TR 100] 2.36x104 = 2.13x10" = 164x104 = = 6. 19x10? ~=— -1..92x103 
200] 1.11x10 9.51x10*  6.58x104  1.22x10% 5 .08x10* 
a/h ‘(iaws $$ —__—— 
2 = ——_— ———- ——___— aft< 1 
10] 1.50x2 1.42x10* 1.05x10°  §.63x10* ~=—-:1.07x107 i i ; sae 
20! 5.27x102 4.99102 «3.73x202S «180x107 Ss: 3.602. lo] 6.15x10> 6.0310 5. 98x10 5.3101 3.08x10" 
50| 2272x102 2.58102 =-1.87x102 = 8. 46x10" Ss. 84x10 ” 2.94x10, piginragea 2.81110" rapioaes fj oma 
100 6x10? -«-8.80x102 = 6.28x102—S «281x102 ~— 6.17107 = Se 2.27R10! —2.22x10?2.04x10° appa! 
200] _ 313x104 x10 2, 1204 bx10? «1. 80x 0* loo} =1.11x10* + 1,0axi0o* 1.04x10 9,10x103 = -2. 05x10 
sts FRE once Eee ld SO occa 200] 5.27x10  5,12x10 4e81x10 3691x104 = -1.67x104 
a2 500] 4.09x10?  3.92x10? 3. 52x0° 46x10? =. 38x10 
E goo] 1.17x10° 1.12x10 = 9.60x109 ~=- 5.96109 -2. 37103 
1200] 2.89x10® + .2,69x10° §=—-2.26x10® 9 1.23110 9-6, 52x09 
1500} 474x10° = 4. 38x10 = 3.59x10° = -1..78x10° ~=s_ +8. 79x10 
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1/40) 


x 
| 
TABLE 
ranks Variation of og,/(P*/a?) with (a/h) and ¢ for x 0 s/L = 
1/40 
x 7 pee ~ — : —~ ae -—— TABLE 12 
° ?, > 9 ; , 
SS Be os Ee. ms a Variation of % Fa*) with (a/h) and x for ¢ | 
ss Fe atte ° 1/40 
—_ a/L = 1/5 
5617x1072 a : > es - ; 
eo 1.50x10“ 9 .02x10™ 5.45x10~ 1.71x10~ -1.1ox10° x x] x2 x3 x) 
: a .27x107 291m? 52x20? Ss .59xl0h = 2, 21x10* SS eS 
223 a aot 
bdo" .72x102 «1627x102 ss 4. 63x10“ +. -1.97x10* = -2.. 30x10+ a/h ee 
bal 3 a : 
08x02 .26x10° 256x102 6. 75x10" 9 -1.14x103 = 8, 50x10* ; 
eJOXLV if 3 2 3 S| fy "7 88&x10* 
13x10 9.19x10? =5.45x10* -4.22x102 = -2.69x107 1 88x10" 
-——_—-— — —— 20 0 4.93x10* 
a/L=1l 50 | 5.11x103 
,81x10~ 2 2 2 1 100 2.75x10 
663x102 2.27x10 1.67x10 1.30x10 8.86x10 -2.80x1072 200 1. 778105 
839x102 2.46x10  6.06x102 -4.60x10* 2, 95x10? - 3. 96x10" 
128x102 4675x102 3.26x103 2.34x10? 1. 33x10? -4.45x10° a/L = l 
96x103 1.73x10  1.13x10  =7.70x103 ~=-3.72x02-S ss -1. 81x10? - 
: / 2 0 
6.210% 3,83x10 = 2. 40x10 = 8.72x107 =5..47x107 we 
30x10? «1482x107 9.52x104 = 1.,.02x104 = -1. 38x10 » : 
67x1 3.91x10 1.77x10° -1.88x104  -1.71x104 100 
58x10 —-7.42x10? = -2. 78x10? = -1.07x10 = - 1. 511x104 500 D 
35x10°  —1.05x10° 3.40x10°  -2.06x10° -8.20x10? 1500 : 
(6/L = 
TABLE 10 PaBLe 13 
4 Variation of o,/(P*/a with (a/h) and ¢ for x = 0; (6 L = Variation of ¢ V,* with (a/h) and x for ¢ = 0 és/L = 
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curves which are given in Figs. 3 to 6 approximate 
lines [particularly in the lower region of 
indicating that simple empirical formulas re- 


straight 
(ah) |, 
lating each of the quantities to a power of (a h) could 
be obtained for the fixed values of (a L). 

Although all stress and displacement quantities thus 
increase monotonically with increasing (a /), no such 
general observation can be made with respect to changes 
in (a L). However, the data do show that within 
the ranges of (a L) and (a h) considered, the flexi- 
bility w (P* Ea’) (aL) 
for fixed values of (a h), while the flexibilities w,, 
(M,* Ea*) and w,, (.\/,* Ea*) as well as the bending 
stresses increase as (a L) increases. For (ah) < 200, 
the radial flexibility is influenced much more strongly 
by changes in (a L) than are the other flexibilities. 
The results also show that the bending stresses as well 
as the axial and circumferential flexibilities are rela- 


decreases with increasing 


tively insensitive to changes in (a L) when (a L) goes 
from 5 to 10, particularly in the region 10 < (a h) < 
200. 

Figs. 3 to 6 also show that the circumferential bend- 
ing stress is never less than the axial bending stress, 
and that the ratio of the former to the latter increases 
as (a ht) increases. On the other hand, the circumfer- 
ential membrane stresses are larger than the axial 
membrane stresses for all (a #) when (a L) = 5 and 
10, and less than the axial membrane stresses at the 
upper end of the (a /) scale when (a L) = | 5 and 1. 
Furthermore, the membrane stresses are generally con- 
siderably smaller than the bending stresses, but never- 


theless, they are not insignificant. 
Displacements and Stresses Along Coordinate Axes 


The quantities discussed above represent maximum 
values with respect to location in the cylinder, since 
they correspond to positions under the load resultants. 
Additional information on the behavior of cylindrical 
shells under localized loadings can be obtained from 
Tables 2 to 21. In these Tables, data for pertinent 
stresses and displacements are presented correspond- 
ing to distances along the axial and circumferential 
coordinate axes equal to one-half, one, and twice the 
half load-length 6 (see Fig. 2) and to one-half the dis- 
tance between the location of the load resultants and 
The corresponding non- 
dimensional coordinates are tabulated in Table 2 for 
values of (a L) equal to | For (a/L) = 
1 5, data are presented for each of the three loading 
cases for (a h) = 10, 20, 50, 100, and 200; for (a' L) = 


one end of the shell (Z 4). 
5 and |. 
|, similar results are presented corresponding to (ah) 


= 20, 50, 100, 500, and 1,500. It should be noted 
that in view of the assumed simple supports at the 


cylinder ends, all radial displacements as well as nor-- 


mal and membrane stresses vanish at x = (L/2a) for 


all loading cases. 
Application of Results to Cylinders Under Arbitrary Loads 


The practical utilization of the data presented in the 
Tables or curves is neither limited by the specific line 
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Fic. 7. Methods for approximating line and surface loadings 


loadings analyzed nor by the fact that in the calcula 
as eZ IQ). For 


example, if the load resultant in a given case is not at 


tions the half load-length 6 was fixed 


the center line of a generator, the present data can still 
be used because of the localized nature of the results 
This is particularly true of the bending stresses, and 
these are usually the most significant quantities of in 
terest. A particular example worked out in reference 
| for the P* case showed that for a fixed 6 L 40 the 
maximum bending stress decreased by only 7 per cent 
when the load was moved from the midpoint of a 
generator to a position equal to SO per cent of the 
distance between this point and the end of the cylinder. 
With these observations in mind, it is not difficult to 
develop a procedure for applying the present numerical 
results to the approximate analysis of cylinders in cases 
for which the actual load-length is not equal to that 
considered herein. 

If for P* or \/,* loading cases [Figs. 2(a) and 2(c) | the 
actual load-length 26* 
need merely be subdivided symmetrically about the 


is greater than 26 = L 20, it 


origin into an odd number of parts each of length L 20 
lsee Fig. 7(a)]. If, after this is done, remainders exist 
at the ends of the actual load-length they can either be 
neglected or approximated by lengths L 20, depend 
ing upon which procedure produces the better ap- 
proximation to the actual load-length, or they can be 
accounted for with the use of a procedure described 
later. Thus, with the loading so subdivided, each 
seginent of length 26 = L 20 can be considered sep- 
arately as though it acted symmetrically at the center 
of the cylinder, and hence, the results given in Figs. 
3 to 6 and Tables 2 to 11 or 17 to 21 can be used. Of 
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TABLE 20 
Variation of —o W.*/a*) with (a/h) and ¢ for x 0; (6/L = 
1/40 
0 g 4 oe ? 
Ng 0 1 Po ?3 Py, 
I em 
a/h a/L = 1/5 
1¢ 0 1.47x102 2.26102 3S 2. 7008 3429 
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é 
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r r >> 5 “as 2 
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TABLE 22 
Approximate Calculations Corresponding to a Uniformly Dis 


tributed Radial Load Over a Square Surface Area 26 X 26 and to 
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Comparison of Present Calculations With Those of Reference 11 
for a Uniformly Distributed Radial Load Over a Square Surface 
Area 26 X 26 
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course appropriate recognition must be given to th, 


symmetric or antisymmetric properties of the dat 
0 planes. Thy 
ther 


with respect to the x = 0 and ¢ 


“unit” 


data obtained for a typical loading cai 


be plotted or tabulated in accordance with its actual] 


relative position along the x axis. When this is don 
for each of the unit loadings, the principle of super 
position can be applied to determine the results of th 
complete loading. Because tabulated data in 
are given only for x 0, for a thorough analysis, 
procedure requires the estimation of the stresses in th 
vy¢ surface from the data given along the coordinat 
axes. 

The procedure described is not limited to uniform] 
distributed line loadings, since by choosing the loa 
intensity of each unit loading in an appropriate manner 
any arbitrary line loading can be approximated. 
thermore, in accordance with previous observations, the 


loading need not be located at the midpoint of a gener 


ator. Thus, general loading cases of the antisymmetri 
\,* type [Fig. 2(b)] are amenable to approximat 
analysis by this procedure. However, here it would 


seem sensible to make use of the tabulated results of the 
JI.* case (Tables 12 to 16 
tributions corresponding to the half unit load-lengths 
The 


remaining portions of the actual load length would 


in determining the con- 
on the positive and negative sides of the x origin 


then be approximated by an even number of unit 

lengths (1 20), and hence, the results of Tables 2 to 1] 
see Fig. 7(b 

26 L/20 


Fig. 7(c) ], the above procedure must be modified. 


could be applied as described previously 

If the actual load-length is less than 
The 
actual loading can be replaced by symmetrical series 
of an odd number of positive and an even number of 
negative unit loadings, with the set of positive loads 
starting at the origin and the negative loads starting at 
distances equal to one-half the actual load-length |see 
Fig. 7(d) J. 
to the actual loading plus two symmetrical loadings 


In this manner the net loading is equivalent 


over lengths equal to one-half the difference between 
the unit and actual load-lengths [see Fig. (7e)]. These 
latter loadings can be made to appear at distances 
from the origin sufficient to eliminate their effect near 
the origin by choosing a suitable number of positive 
and negative unit loadings. In fact they may be trans- 
ferred to positions quite close to the end supports of 
the cylinder. This procedure can be combined with 
the previous one for all load cases and may be par- 
ticularly useful if it is necessary to choose small inter- 
vals when approximating variable load distributions. 
The present results for line loadings along generators 
can also be applied to loadings distributed in the 
as axial directions (surface 


circumferential as well 


loads). The given loading must be subdivided into 
longitudinal slices, each of small finite circumferential 
width. The loading represented by each slice can 
then be approximated by a statically equivalent line 
load whose profile in the x direction approximates that 
The effects of this equivalent load can 
After 


of the slice. 
then be determined as described previously. 
each slice has been considered separately, the data ob- 
stresses—can be plotted or tabulated rela- 


tained—e.g., 


planes 
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Dynamic Considerations Relating to the 


Behavior of Inertial Space- 


Stabilized Platforms 


EDMUND 


GORCZYCKI* 


Sperry Gyroscope Company 


INTRODUCTION 


oo THE HUMAN ELEMENT is present to 
close the multiloops in a moving vehicle such 
as an airplane, corrections are provided manually to 
steer the craft to its destination. This is especially 
true if the airplane is flying at subsonic speeds, since a 
navigator can determine by various means his present 
position and his course, taking into account natural 
effects such as winds and gusts. He has sufficient time 
to make computations to ascertain the course at any 
time to arrive at the touch-down point. 

With the advent of unmanned aircraft, such as guided 
missiles, or manned aircraft at supersonic speeds, a 
requirement exists to derive basic information regard- 
ing aircraft altitude, rate of change of altitude, rate of 
travel at any time with respect to the earth, and many 
other time variables which a navigator has to compute 
in a subsonic aircraft. To meet this requirement, 
the inertial space-stabilized platform was conceived. 
Such a platform is supported within a housing which 
is rigidly attached to the body structure of the air 
frame. 

The subsequent discussion presents the dynamics of 
an inertial space-stabilized platform in the presence of 
an arbitrary housing disturbance due to air-frame mo- 
tion at an isolated point in inertial space. Air-frame 
translational motion with respect to the earth is not 
considered. 

This paper avoids duplication of material regarding 
explanation of gyro operation in detail since it is very 
well described elsewhere in the literaturet. The stress 
in this paper lies in the three-dimensional aspect of 
platform operation to determine the coupling effects 
which influence the platform motion. The basic fac- 
tors which must be considered in stable platform de 


sign are discussed. 
GENERAI 


The inertial space stable platform is a device which 
has the capability of maintaining its vertical aligned 
with the direction of local gravity independent of air- 
frame rotational and translational maneuvers. To 
satisfy such a requirement, the device must behave as 
an integrator about each of three mutually perpendicu- 


Presented at the Instruments Session, Twenty-Fourth Annual 
Meeting, IAS, New York, January 23-26, 1956. 

* Senior Engineer. 

t See footnote, p. 138. 


lar reference axes in inertial space. The reference axes 
are established by three integrating rate gyros whos 
input axes are mutually perpendicular. The input axis | 
of one of the gyros is aligned perpendicular to the } 
platform while the input axes of the remaining tw 
gyros are located in the plane of the platform and arbj 
trarily oriented. 

The gyros serve the function of detecting the platfon 
disturbance and provide a reorientation of the inertial 
reference frame in accordance with prescribed command } 
signal data. The output of the gyro is a rotation about 
the gyro output axis proportional to the time integral 
of the angular velocity about the gyro input axis 
Application of a torque about the gyro output axis cor 
responds to an equivalent rate about the gyro input | 
axis. This action results in a rotation of the input axis } 
about the output axis to effect a change in reference | 
coordinate. 

Stabilization of the inertial space platform is accom | 
plished by closing servo loops around the individual | 
integrating rate gyros. Servo drive motors rotate the 
platform and its associated supporting gimbals in re 
sponse to the electrical signals derived from the gyri 
rotation about its output axis. The motor drive is 
continued until the output signal from each gyro is 
nulled. When this null condition is satisfied, the plat- 
form is returned to the same space orientation that it 
had prior to the disturbance which upset it and the 
amount of platform disturbance is indicated on an elec- 
trical pickoff device mounted at each of the gimbal axes 


and the azimuth axis. 


ANALYTICAL DEVELOPMENT 


The Newtonian law of motion for a rigid body is ex- 
pressed generally in concise vector notation as follows 


MP fs = m[Pins +2 Xtm taxr+ax 
(a xXFf)] ( 
where 
7 = any directed distance referred to either 


fixed or moving space 
angular velocity vector representing rota- 
tion of the rigid body about an axis ol 
rotation passing through the body 
= behavior of the distance 
respect to moving space 


(  )ms directed with 


130 


pone 


ire I 


\PPI 


towd 
side 
form 
ixeS 
are 
stab] 
coor’ 
cent 
the | 
initic 
are Z 
XxX} 
¥ 
Z; 
The | 
TI 
its a 
axis | 
refer 
troid 
locity 
and 
along 
Th 
the | 


consi 


e 


rence axes 
ros whos 
input axis | 
ar to th 


ining tw 
and arbj 


platforn 
ie inertial | 
command j 
ion about 
e integral 
put axis 
AXIS Cor- 
ro input 
nput axis 
reference 


iS accom 
idividual 
otate the 
ils in re- 





the gyro 
drive is | 
| gyro 1s 
the plat- | 
n that it 
and the 
| an elec- 
nbal axes 


dy 1s ex- 
follows: 


o either 


ng rota- 


) axis Ol 
dy 
ce with 


corresponding behavior of the directed dis 
tance with respect to fixed space. 


If a right-handed reference coordinate system is 
chosen as shown in the insert of Fig. 1, the components 
of acceleration of any moving point 7? in the platform 
with respect to fixed reference axes, Y;, V7, 27, are ob 


tained from an expansion of Eq. (1) along these axes 


is follows: 
oe 
'20, — VW. — "ry\W 7 @ i 
Vr WyrW, TT XeyB (2 
} + 2lw.(7 —- w,(F “a 
VW - Fie — FAG +r W + 
TrWrW, TT yw (4 
} 7 2[w,(7 — W\Tr)r | — 
Vr ,@ = by — TAG TT Wy a 
} WW T VyjWye (4 
Eqs. (2), (3), and (4) shall be used in the develop 


ment of the force equations acting on a point mass m™ 


in the platform. Suksequently, these force com 
ponents are used to derive the torque equations which 


are required in determining the motion of the platform. 


\PPLICATION OF NEWTONIAN FORCE EQUATIONS TO A 
PITCH-ON-THE-OUTSIDE STABLE PLATFORM WITH 
STATIONARY HOUSING 


Without sacrificing generality of exposition and 
mathematical development, regarding stable platforms, 
in general, the subsequent discussion will be directed 
toward the dynamic behavior of a pitch-on-the-out- 
side stable platform whose configuration in rudimentary 
form is shown in Fig. 1. V5, 245 


axes in the insert of Fig. | are fixed in inertial space and 


The reference X,, 


are aligned initially with corresponding axes in the 
stable platform. It is assumed that the origin of the 
coordinate system lies along the azimuth axis at the 
centroid of the platform. The positive directions of 
the platform reference axes are the same as those of the 
initial airframe reference axes, when the gimbal angles 


are zero, that is, 


X; positive toward the nose 
V; positive toward the right wing 
Zi positive toward the bottom of the air frame 


The subscript / designates reference to inertial space. 

The platform is assumed to be oriented initially with 
its azimuth axis aligned with the fixed Z; reference 
axis of Fig. 1 and the plane containing the X; and ¥, 
reference axes passes perpendicularly through the cen- 
troid of the aximuth axis. The platform angular ve- 
locity of magnitude G;_, is assumed not to be constant 
and the platform angular velocity vector is directed 
along the positive azimuth axis direction. 

The accelerations of a typical particle P of mass m in 
the platform located at coordinates (X, Y, Z) shall be 


considered for the following conditions: Rotation of 
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the platform about its azimuth axis is the only motion 
involved. The azimuth axis is assumed to lie along the 
reference Z,; axis and rotation is in the positive sense 
as shown in Fig. 1. Then, in addition to the rotation of 
the platform about its azimuth axis (Z; axis), assume 
that the platform has an angular velocity G;_;¢ about 
the inner gimbal axis (platform and housing Y axis 
initially) and an angular velocity G;_o, about the outer 
gimbal axis (platform and housing } axis) in the indi 


cated positive directions of Fig. | 


The composite accelerations at point ? (x, y, 2) due 
to simultaneous rotations about the Y,, }7;, and Z; 
axes are derived from application of Eqs. (2), (3), and 
(4) and are given by 
(7, )y —(G*, p 7 G";_o0Gg)X G; OG 

Gi rG WG oaN ) 
CAUP —(G*; G*; iG}y > Gi_1G2 4 
Gy co G& GG] oOGN 10) 
tRP (2G, Gr _ Gr OG)X tT 
(Gi_1g + 2G, oGgGl p)¥ 
(G*, 1G 7 G*; OG)Z (7 
Fig. 2 is drawn on the basis of Eqs. (5), (6), and (7 


Since point ? is only a single point out of the multi 
tude which comprise the rigid body (the platform 
using Eqs. (5) and (7), and a double subscript notation 
for principal moments of inertia and products of in 
ertia, the effective couple (7%) about axis OY, is 


T 0G LG og — UI G, » Gi 19 = as G1_1G 
(21, )G) OG IG, — | G 
1,(G?, _— G?, IG T l, Ge oh OG (S 


Similarly, using Eqs. (6) and (7), the effective couple 


(Tyq) about axis OX, is 


T 1 T-:G1 IG = (/ tne Gs OG = PE er OG T 
((27,,)Gr 1G Moe a I, Gr pT 
1,(G%1-p — G*1~-0g) — 1,:G1-1cGr-oc (9 


Using Eqs. (5) and (6), the effective couple (7) about 
axis OZ; due to motions of the platform about the inner 


and outer gimbal axes is: 


ie G1 >= [,Gi-1qG — [,:G1-0G — 


/ 


I.,(G*: io G*; OG (10) 


Eqs. (8), (9), and (10) are represented diagrammati 


cally in Fig. 3. It is seen from these equations that in a 
dynamically unbalanced stable platform, motions of the 
platform about the inner and outer gimbal axes mani 
fest themselves as an effective torque about the plat- 
form azimuth axis. Similar dynamic results are ob- 
tained with respect to the inner and outer gimbal axes 
due to angular rotations about the azimuth axis. 

If the stable platform housing is held stationary with 
respect to the established stationary set of reference 
axes, the outer gimbal axis is fixed in the J}; direction. 
The aximuth axis of the stable platform has two de- 
grees of freedom with respect to its supporting gimbals 








due to its construction. Thus, if the outer gimbal 
axis is displaced by some angle G;_;¢ with reference 
to the Z,; axis, the outer gimbal angular velocity 
G;_og Will no longer be perpendicular to the azimuth 
axis. However, the entire preceding mathematical 
development is based on the fact that the angular ve- 
locity vectors lie along three mutually perpendicular 
fixed space axes in the plane of the platform (but 


not rotating with the platform Thus, to validate the 


mathematics to encompass other situations, it is 
necessary to resolve the angular velocity G;_oq and 
angular acceleration G;_oq vectors into two components. 
That is, one component of each of these vectors will 
be parallel to the azimuth axis and the other component 
perpendicular to the azimuth axis. The components 
Gene 0 Ge.se nl Gece Ge Greco parallel to the 
azimuth axis play no activé role in the determination 
of platform motion about the platform inner and outer 
gimbal axes and do not affect azimuth motion. Thus, 
the effective components of outer gimbal angular ve- 
locity and angular acceleration are G og cos G;_7g and 
Grog cos Gy_1¢, respectively. Combined with the 
effects of gimbal inertias, these facts are used to derive 
Eqs. (15), (14), and (15), which represent torque equa- 
tions relating to a pitch-on-the-outside stable platform 
with a stationary housing angularly displaced by an 
angle G;_;¢ with respect to the stable platform. 

The development, thus far, did not include the effect 
of the moments of inertia of the gimbals which support 
the platform. By reason of the gimbal construction 
and housing orientation within the air frame, the inner 
gimbal axis is always in the plane of the platform, re- 
gardless of the housing orientation in inertial space. 
Defining Jjgr as the moment of inertia of the inner 
gimbal frame about the inner gimbal axis plus asso- 
ciated gearing, electromechanical components, static 
balancing weights, and any other masses associeted 
with the inner gimbal frame, the torque required to 
accelerate the inner gimbal frame by itself about the 


inner gimbal axis is expressed by 
1G lie: GI 1G 11) 


Define Jogr as the moment of inertia of the outer 
gimbal frame about the outer gimbal axis plus associated 
gearing, electromechanical components, static balancing 
weights, and any other masses associated with the 
outer gimbal frame. 

The moment of inertia of the structure (exclusive 
of the platform) which it is desired to accelerate about 
the outer gimbal axis is not only Jog but also some func- 
tion of gimbal displacement between inner and outer 
gimbals acting on /;gr. Assuming that the platform 
housing is oriented such that the inner and outer gim- 
bals are perpendicular to each other, the moment of in- 
ertia of the combined support structure mass about the 
outer gimbal axis is Joge + Jjgr. In general, the inner 
gimbal frame will be a shell-type ring. The moment 
of inertia of the inner gimbal ring about the outer 
gimbal axis is reduced from J;ge when G;_7g = O° to 
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Tigr/2 when G;_1¢ G0) Based on this discussjo; 
the torque required to accelerate the combined inner 
and outer gimbal frames by themselves about th 
outer gimbal axis is expressed by: 


’ ‘ ( . 1 
Lies COS Gy Ie 2 { Gi ( iZ 


Too = { Toa: 


Eqs. (8), (9), and (10) can now be modified to i 
clude the effects of gimbal inertias and housing tilt 


follows: 


: ‘ ex 
Toc ’ Tog: T lias cos Gry_7¢ 2 ; 


cos Gy; 1G) COs Gy; IG 


s 
| 
+ 
— 
~ 
~ 


cos Gy, 


G?, 1g) cos G; re ai | [Gy 1GGr Or cos” Gy; 


\ 7 res + 1 Gy] lt / 1Gr og COS G) 
' 1G nee B (J Gy; Gy og COS Gr I 
2] y Gy; 1G IG] T G I 9a G*; og COS 


/ G) GQ OG cos G, I¢ 


7, / Gy I = ,. Gi IG / GI On cos G; Ie 


G*; 1 — G*;~0G Ccos* Gr_ 1 [ 


Neglecting frictional effects, Eqs. (13), (14), and 
are the differential equations of motion of the pitch-on 
the-outside stable platform (including the effect of | 


supporting gimbals) with a stationary housing. Simul 
taneous solution of these equations with proper initi: 
conditions yields the results of platform rotation Gy_¢ | 
about the outer gimbal axis, G;_7;¢ about the inner 
gimbal axis, and G;_, about the azimuth axis. 
From the dynamics of a rigid body in inertial space, | 
the rate of change of angular momentum about a! i 
axis of rotation of the system is equal to the torque ans nia 
ing from the impressed forces about the same axis 
Thus, with respect to the torque equations (13), (J4 
and (15), each term on the right-hand side represents 
a rate of change of angular momentum. The bol 
5, and 6 show the angular 


linear vectors in Figs. 4, 
momenta involved in arriving at the derived torques 
of Eqs. (13), (14), and (15) about the Y;, ¥7, Z; axes 

respectively, when the time rate of change of thes 

angular momenta is taken for a platform with stattot 

ary housing. The thin-lined vectors in Figs. 4 and °| 
represent the components of the bold-lined vect rs | 
which are projected along the X,;, };, and Z; axe 
assuming that the platform vertical is disturbed fro! 

the Z; axis during platform stabilization by smal 

angles «G;_og and e«G;_;q about the outer and inner 
gimbal axes, respectively. 

The angular momenta of the gimbal frames abot! 
the X; and ¥; axes, exclusive of the platform inertia, at 
Teoshes . rg and (Tor + Jigpr cos? (Gy IG 2) [cos Gy 
Geuae respectively, and are directed in the same sens aa: 


as the respective corresponding gimbal torques. 
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PLATFORM TORQUE EQUATIONS DUE TO 


ROTATING HOUSING 


EFFECT ON 


Thus far, the development considered only the ef 
fective torques acting on the platform whose housing 
is stationary. It is now desirable to investigate the 
effect on the platform torque Eqs. (13), (14), and (15) 
contributed by rotation of the platform housing in 
inertial space. From observation, it is apparent that 
the only rotation of the platform housing which will 
effect a physical rotation of reference axes fixed to the 
platform gimbals and aligned initially with the Y;, 17 
axes is a rotation about the inertial Z; axis. 

The rate of rotation of the housing about the Z; axis 
in inertial space shall be designated as Gz, , for which 
the quantitative expression 


Gz, =n —psin Gy,_, sec Gy; 4 + 


r cos Gy,_,, see Gx,;—, (16) 
may be derived. 

The impressed torques Zo, 7Zic, and 7), of Eqs. 
(13), (14), and (15), respectively, give rise to the angu- 
lar velocities G;—o¢, Gr—1q and G;_, about the plat- 


form outer, inner and azimuth axes when the housing 
is stationary. Acting upon the inertia constants of the 
platform, these angular velocities give rise to angular 
momenta along the Y,;, ¥7, and Z; axes. 

In the process of platform stabilization against physi- 
cal disturbances, the torques impressed on the plat- 
form act in such a direction as to counter the disturbance 
about the axis which is affected by the disturbance. 
If the stabilization of the vertical is not exactly perfect, 
the platform vertical will be displaced with respect to 
the and 
€G] 
As a result, components of angular momenta of the 


true some small angles eG; —0¢ 


1 about the outer and inner gimbals, respectively. 


vertical by 


platform in Fig. 6 are transferred so that they align 
with the original angular momenta under consider- 
ation. 
fixed to the platform gimbals are rotated at the rate 
G.,—-» and cause the angular momenta along these axes 
The rotation of these 


Through physical constraint, the reference axes 


to be rotated at the same rate. 
angular momenta gives rise to a rate of change of 
angular momenta (or torques) in the same sense as 
the angular rotation. Similarly, angular momenta of 
the platform in its gimbals contribute their rate of 
change to affect the torques which act on the platform. 
Fig. 7 shows these additional torques on the platform 
which arise from platform housing rotation. From 
Fig. 7, the following equations can be written directly. 


AT 04 [/ G; v| leGr—oa| [cos Gr_1¢] [Gz, al + 
WrerGr 16) [Gz,-H| cos Gr—1g — 
(, G1 —16 | [eGr—oa | [Gz, uw} cos Gr_1e — 
(1,.G; al (Gz, w| cos G;~1G — 
(2 -yGr og cos G;_—14] (Gz, nw| cos Gr_1g — 


[Gr og COS Gy _16¢ I[eGy og |[Gz, u| cos G;_ 16 


AGRONAUTICAL 
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AT 16 = [[.4G1—p] [eGr—1¢)[Gz,-4) + WyGr—-p] (G1, 9] 4 
[21 ,,Gr 14 | [Gz, al — [, Gr—o« cos Gy 14G| X 
[eGr—1¢] [Gz,-n| — Us:Gr—1e) [eGr—16] X 
Gz, 7 [Loar t Lie cos? (Gy —1¢/2 | x 
[Grog cos Gr—16| [Gz,—» Is 
AT, 0 19 


Eqs. (17) and (18) represent the additional torqu 
components accrued along the outer and inner gimbal 
axes, respectively, due to platform housing rotatior 
about the inertial Z; axis at angular rate Gz, wu. Eg 
(19) indicates that no change in torque takes plac 
about the platform azimuth axis due to platform hous 
ing rotation. Eqs. (17), (1S), and (19) must be added 
algebraically to Eqs. (13), (14), and (15), respectively, 
to obtain the total dynamic behavior of the platform 
(including the effects of gimbal inertias) in the case of 
the rotating platform housing. 

Fig. S is a block diagram representation of this sum- 
mation of torques, showing the input torques at the 
platform about the platform outer, inner, and azimuth 
axes and the output platform angular accelerations and 
velocities. The expressions for the torques experienced 
at the platform are entirely general and any degree of 
simplification can be made upon them, depending upon 
platform design. For instance, by elimination of ap- 
propriate terms the platform may be transformed from 
a dynamically unsymmetrical one to a dynamically 
symmetrical one. Inclusion of the A7’ terms permits 
study of the effect of rotating housing compared to the 
case in which the housing is stationary. 

The following development provides a transforma 
air-frame axes to 


tion of rotational velocities about 


equivalent rotations about inertial space axes. From 

this derivation, an expression is derived for the rota 

tional rate Gz, a Which is required in Fig. 8. 

TRANSFORMATION OF ROTATIONAL VELOCITIES ABOUT 

AIR-FRAME AXES TO EQUIVALENT ROTATIONS ABOUT 
INERTIAL SPACE AXES 


Since the stable platform housing is rigidly attached 
to the air frame, its rotational velocity in inertial space 
is identical to the rotational velocity of the air frame 
Air-frame angular velocities, and hence the platform 
housing, about air-frame body axes may be trans- 
formed into equivalent rotational velocities Gx, Hy 
Gr, y, and —_— about the inertial space axes X), 
}7, and Z;, respectively. 

Air-frame angular velocity relationships based on the 
sequence of rotations yaw, roll, and pitch may be 
developed by matrix methods or by conventional 
methods of resolving rotational velocities represented 
by vectors. The magnitude of the rotational velocity 
vector is determined by the amount of rotational ve 
locity, and its direction is obtained by means of the 
right-hand rule. 

Referring to Fig. 9, the original orientation of the 
air-frame axes is represented by a mutually orthogonal 
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Fic. 8. Block diagram of gyro mounted stable platform. It is assumed that the input axis of the roll integrating rate gyro is initially 
aligned with the inner gimbal axis, and that the platform is initially erect 
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frame of axes with ‘‘O” subscript, 1.e.—.Yo, ¥o and Z). 
The first sequence of rotation is a yaw maneuver at an 
angular velocity 7 about the Zp) axis. This action dis 
places the air-frame axes from their original orientation 
into another axes orientation designated by the ‘] 

subscript i.e., i, i, Z;. At this time, the subscript 
‘1’ frame is rotated at an angular rate £ about the Y, 
The roll maneuver rotates the }j, Z; axes into 
while axis Yo is still aligned with axis 


AXIS. 
position V2, Zo 
X,. The final maneuver rotates the axes with subscript 
2" into position with subscript ‘‘3"’ at an angular ve 
locity f about the FY» axis. 

On the basis of Fig. 9 and the preceding discussion re 


garding air-frame motion, it may be shown that 
P é Ccoo« — % sin ¢ cos & 20) 
q €+ nsiné 2] 


y €sin ¢ + n cos ¢ cosé (22 


where p, g, and r are rotational rates about the air- 
frame roll, pitch, and yaw axes, respectively. 

By physical reasoning, it can be deduced that for a 
pitch-on-the-outside stable platform, the sequence of 
rotations through which the platform housing must 
be rotated in order that each of the three angles about 
the three mutually perpendicular platform axes are 
independent of each other is yaw, roll, and_ pitch. 
Thus, the platform housing orientation is exactly the 


same as the air-frame orientation. Therefore, 


. . (92 
5 Gxj-4 \<9 
¢ Gy, F (24 

, ar 
Ui] Gz, H 49 


Substituting Eqs. (23), (24), and (25) into Eqs. (21 


(21), and (22) and solving, 
go F sin Gy, H (26 


Gx,—H p cos Gy, 


Gy, H Pp tan Gx, H sin Gy — + G=— 
r tan Gx, H COS Gy, H (27 
Gz; H =—— sin Gy, H sec Gy a i 


r COS Gy H See Gx! H 28) 


Eqs. (26), (27), and (28) are the required transforma- 
tions of platform housing angular rates from air-frame 
body axes into equivalent rotational velocities about 
the inertial space axes X;, ¥7, and Z;, respectively. 


DETERMINATION OF AIR-FRAME INITIAL ANGLE 
CONDITIONS 


In order to specify a unique orientation of the air 
frame in space, it is necessary to relate the air-frame 
axes to reference axes. This is accomplished by means 
of a set of three air-frame orientation angles (y, 6, o), 
called Euler angles, and a choice of the sequence of 
rotations which the air frame must follow by the 
amount of Euler angle corresponding to the air-frame 
axis with which it is associated. The choice of rota- 
tions recommended by the NACA is yaw y, pitch 6, and 


ICAL SCIE 


NCES F 


roll @, as shown in Fig. 10. 


had been chosen in what follows so that the amount of 


Euler pitch angle in space with respect to the referenc 


coordinate system in the Y¥; — ¥, plane is unaffect 


by either the Euler yaw angle or Euler roll angle 


By equating elements of the air-frame orientation 


matrix with corresponding elements of the platfor 


Orientation matrix, the following results are derived 


Gx; H sin cos 6) sin o mA) 
Gy H tan tan 6) sec @ HU 
Gz,-H Yo — tan tan @» sin 6 3] 


where the subscript “‘O’’ notation indicates air-fram 


and platform housing initial orientation angles with 


respect to the reference inertial axes Y,;, 7, Zr. 


The advantage of deriving Eqs. (29), (30), and (3 


in the sequence of rotations yaw, pitch, and roll is that 


this is the sequence generally applied in aerodynami 
analyses. 
If a perfectly stabilized and erect pitch-on-the-out- 


side platform is considered, with no friction at the 


pivots, Eq. (29) represents the indicated angle between 
the inner and outer gimbals, Eq. (30 
indicated angle between the outer gimbal and housing, 
and Eq. (31) represents the angular rotation of the 
housing about the platform azimuth axis. 


STABLE PLATFORM ISOLATION CONSIDERATION 


Eqs. (13), (14), (15), (17), (18), and (19) show the 
relationships between the torques acting at the plat 
form and its gimbals and the resultant motion that 
evolves from the presence of these torques. The dri 
ing torques originate in the motors at the gimbal and 
azimuth axes. 

The motion of the platform gimbals with respect to 
inertial space can be related to the servo motor shaft 
torque, housing motion, interfering moments consti- 
tuted by stiction and external loading. If a free body 
diagram is drawn for all the gears between the servo 
motor to the platform, it can be shown that 


Ttoqg + (f (gear inertias, motor inertia, vis- 


Tog = 


cous motor damping ratio)| (Gy,-4) —  [g(gear 


inertias and ratios, motor inertia, viscous motor damp 


) 


gear interfering moments (32 


ing (Gi OG) — 2 


Equations of similar form can be derived for 77¢ and 
i> 

The equations for 7Zoc, 77¢, and 7, in the form of 
Eq. (32) constitute dynamic transformations of the 
torques at the servo motors which drive the platform 
gimbals to maintain the verticality of the platform and 
to orient a set of reference axes in the plane of the 
platform. 

Platform isolation, which is the relationship of plat 
form motion to housing motion, ranges anywhere from 
perfect isolation in which no gyro is required for loop 
monitoring to the case in which the gyro-servo loop 
which is used in inertial space stabilization integrates 
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Determination of maximum follow-up rates about the 
pitch gimbal axis 


with respect to the housing independent of housing 


motion. when the 


Perfect isolation is approached 
servo motor inertia, damping, and motor pinion gear 
inertia approach zero; the latter condition is approached 
when the gear ratio from motor to platform approaches 
infinity. 

Normal platform stabilization design is a compromise 


condition between these extreme cases. 


STABLE PLATFORM GIMBAL SERVOS 


It can be shown from the application of servo tech- 
niques to the design of the platform gimbal servos that 
the desired transfer function of platform displacement 
irom the motor shaft at the outer gimbal to the plat- 


form itself is of the form: 


G;_oG\ Pp) 


e(p) = 
K(1 


+ prn)/pUl + pre)(Q1 + prs) (33) 


Similar equations apply to the inner gimbal axis and 
platform azimuth axis. 

Based on preceding discussion, the error signal which 
drives the outer gimbal servo must be altered by the 
secant of the roll angle in order to project it in line with 
the outer gimbal axis. This transformation is shown 
in Fig. 8. 

In order to arrive at a transfer function of the form 
of Eq. (33) the following must be assumed: 

The gimbal 
Plat- 


a) The servo system is linear. (b) 
loading reflected to the motor is constant. (c) 


SPACE 


STABILIZED PLATFORMS 


form angular disturbance occurs about the gimbal 


axis under consideration. Rotations about the other 
two platform axes are zero. (d) The orientation of 
the housing and platform is such that all the indicated 
gimbal angles are zero initially e) The platform hous- 
ing 1s nonrotating. 

In general, the design parameters in Eq. (35) are 
obtained on the basis of relative magnitudes of the con 
stants, practicability of obtaining these magnitudes 
with undue difficulty, and stability and performance 
within the bandwith of frequencies in which the servos 
must respond. The servos must have the capabilits 
of tracking the maximum platform rates determined 
from the prescribed maximum rotational rates of the 
Assuming that the 


about the 


air frame about its body axes. 


maximum rotational rates Pyor, Gmary " 
air-frame body axes and the maximum bank angle 
Omar, Of the air frame are specified, the maximum gim 
bal rates about the inner gimbal axis, outer gimbal 
axis, and platform axis can be determined from Figs. 
11, 12, 12a, 
the sense that in order to arrive at the maximum pitch 


and 13. Figs. 12 and 12a are related in 
gimbal rate it is necessary to divide the ordinate in 
Fig. 12 by the ordinate in Fig. 12a for a corresponding 


value of @, , 
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azimuth axis 


Overdesign of the platform servo may be avoided by 
designing the servo to follow up accurately to some 
percentage of the values derived from Figs. 11, 12, 12a, 
and 13 and permitting the platform vertical to lag 
behind its true position as long as the gyro reference is 
not lost. In a zero velocity error servo, the vertical 
will arrive at its true position in time. 


THE INTEGRATING RATE GYRO 


Expressions have been derived for the dynamic 
motion of the platform with respect to inertial space 
and for the motion of the platform housing with respect 
to inertial space. Both of these motions are “‘open- 
loop” unless a device is placed into the system to derive 
between these two motions. 


the relative difference 


Such a device is the integrating rate gyro. Three of 
these gyros are mounted on the platform to detect 
platform disturbance about their respective input axes 
which are mutually perpendicular. 

The integrating rate gyro is a single degree of free- 
dom gyro immersed in a viscous medium. The gyro 
output angular displacement is the resultant time 
integral of the angular rates of the gyro about its input 
axis. The input axis of the gyro lies along the per- 
pendicular to the gyro spin axis and the output axis. 
The output which is the integral to which reference 
was made is detected using a microsyn pickoff. Fig. 
14 shows a working sketch of sucha gyro unit. Draper, 
Wrigley, and Grohe have presented a comprehensive 
account of the mechanism of the integrating rate 
gyro.* 

On the basis of the dynamics of the gyro, the gyro 


K [pl + Top) |. 
* The Floating Integrating Gyro and Its Application to Geo- 


metrical Stabilization Problems on Moving Bases, Aeronautical 
15, No. 6, pp. 46-62, June, 1956 


transfer function is expressed by 
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Fic. 14 Schematic diagram of the integrating rate gyro 


However, the time constant 7, in (1 + 7,p) term usually 
is adjusted in design so that its effect is small within 
the working range of frequencies. 

In Fig. S, it is assumed that the input axis of the roll 
integrating rate gyro is initially aligned with the inner 
gimbal axis. In general, the input axes of the roll and 
pitch integrating rate gyros are arbitrarily oriented so 
that a resolver is required to transform the platform 
rates in the direction of the gimbal axes as projected 
in the plane of the platform. 


OUTPUTS 


The outputs of a typical system as described in this 
paper are high accuracy indications of air-frame atti- 
tude and rates with respect to an inertial space-stabi 


lized reference. 


CONCLUSIONS 


A mathematical development has been presented 
which describes the behavior of a stable platform which 


is under the influence of air-frame disturbances in 


isolated The method enumerates the factors 


which must be considered in stable platform design. 


space. 


The presentation allows one to develop a deeper in- 
sight into inertial platform operation by considering 
its dynamics; it permits coupling effects between axes 
in the platform to be investigated fully; it establishes 
design goals for the azimuth and gimbal follow-up 
servos based on air-frame performance limits; it shows 
air-frame quantities to 


how the transformation of 


platform quantities is effected. These are a few of 
the many significant benefits derived from an analysis 
of the type which has been presented. 

To extend the development described in this paper, 
it would be desirable to determine the effects of trans- 
lation of the air frame with respect to the earth. Other 
areas of interest lie in determining the effect of tem- 
porary power interruption during an air-frame ma- 
neuverand investigating various methods of establishing 


a platform vertical. 


th 


tid 


OUTPUT 
y AXIS 
(OA) 


e gyro 


1 usually 
1 within 


the roll 
he inner 
roll and 
-nited so 
ylatform 


rojected 


| in this 
ne atti- 


e-stabi 


esented 
1 which 
ices in 
factors 
sign. 

per in- 
idering 
‘nl axes 
blishes 
low-up 
shows 
‘ies to 
few of 


nalysis 


paper, 
trans- 
Other 
f tem- 
ie ma- 


ishing 


The Calculation of Thermoelastic Beam 
Deflections by the Principle of 


Virtual Work’ 


BRUNO A. 


Columbia [ 


SUMMARY 


\n approximate formula for the calculations of the deflections 
f heated beams is studied in this paper. This formula is anal 


the familiar equation v = S (MM>)/ET)dx of strength 


The accuracy of the formula, from the standpoint 


gous to 
f materials 
fa more rigorous thermoelastic solution, is assessed with pal 
rectangular beams under arbitrary tem 


ticular reference to 


perature distributions 


INTRODUCTION 


r I Ve CALCULATION of the deflections of beams is very 
often most conveniently carried out with the aid of 
The deflec 


tion vy at any point in the beam is then usually obtained 


the principle of virtual displacements. 


from the formula 


AM (xo). p(xo, x) 


+] 
v(x - dxo 
a EI 
which is easily derived from that principle. In this 
formula the beam is assumed to extend from x 0) 


to a L; its cross section has a (principal) moment 
of inertia J about an axis perpendicular to the direction 
of v, M7 (x9) is 
the bending moment at any point x = x, and JJ p(x, 


xo by a dummy 


and its material a Young’s modulus F. 


v) is the bending moment caused at x 
system of loads constructed so that it does work if and 
only if the deflection v occurs at x; the magnitude of 
the dummy loads is such that this work is numerically 
equal tov. 

The principle of virtual work is, of course, of general 
validity, and holds regardless of whether the beam de- 
flection is due to heating or to actual applied forces. 
It is therefore desirable to see whether a simple for- 
mula, analogous to Eq. (1), cannot be derived from it 
for the calculation of thermal deflections of beams, and, 
further, whether its range of validity is as wide as that 
of Eq. (1). 
establishment of such a formula. 


The purpose of the present paper is the 


The derivation of an equation of the form of Eq. (1), 
for the deflection of a beam subjected to an arbitrary 
temperature distribution, is an easy matter, provided 
the assumption that sections plane before heating re- 

Received February 24, 1956 
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niversily 


main plane after heating is accepted and that principal 
axes are employed; it takes the form 


[ Ur(x) Mp(xo. x), ~ 
J EI ‘ ~ 


The quantity \/ ris defined by 


. . 
ak] VY, ¥, ZV ay dz Za 


where the integration extends over the entire cross 
prin 


Mr (x 


sectional area A, and where the y and 2 axes are 
cipal axes of the cross section, the y axis extending in 
the same direction as the deflection v. The coefficient 
of thermal expansion is denoted by a. The derivation 
of Eq. (2) must be examined now from a more rigorous 
standpoint so as to determine its accuracy and range ol 
validity. 


(GENERAL FORMULA FOR DEFLECTIONS O! 


DETERMINATE BEAMS 


Let the stresses due to the dummy load system be 
denoted by the subscript D, and let the strains due to 
the actual (thermal) loading be denoted by the sub 
script 7. Then the work done by the dummy system, 
as defined above, will be 
6 SSS (Orp€rT T Fyp€yT 

ToyDYzy7 F 1 
if the deformations are taken to be those caused by the 
actual loadings. The 
chosen so as to cover the entire structure. 


limits of integration must be 
The sym 
bols o and ¢ denote the normal stresses and strains in 
the directions indicated by the subscript, and simi 
larly for the shear stresses and strains 7 and y. 
Consider now a statically determinate beam, under 
an arbitrary temperature distribution and under no 
applied loads; then the reactions will all be zero and 
the thermal stresses will satisfy equilibrium conditions 


of the type 


OrT dy dz | OrTy dy dz 
JA JAe é z 
J.J 


The determination of the total deflection 6 is usually 
most conveniently carried out by first calculating sep 
arately the deflections v and w along the principal direc 


rdy dz QO (4) 
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tions y and g, respectively, of the cross section, and then 
superimposing them according to the formula 6 

to" = Ww") 
will be considered, but it is clear that the analysis pre 
sented can be applied to w with only a change in 
nomenclature. For the present case, then, one may 
write the stress components due to a dummy load 
passing through the shear center and, therefore, causing 


no twist, as follows: 


OrD Mp ay | Tr-D —(O¢ Oy) "i 
T yrD (Oe Oz [| — [y? (2/7 ) |} (d Mp/dxo); ee 
oD o-D TusdD o) 


where the function ¢ must be determined from the 
Saint Venant theory of flexure (reference 2). Sub 


stitution into Eq. (3) gives 


21 _ 
z(x) (1/EI | ila ET — 
0 ie 


V\OyT T O27) |¥ Mp (xo, x) -—(1+ v)yv*(d Alp dx TryT T 


2(1 + v) [(0¢/02) 72,7 — (O¢ Ov) rer] dxy dy dz (4) 
provided that use is made of Eq. (4) and Hooke’s law 


Orr — V(oyr + o:7) + akT tZ 


her 

The first term on the right-hand side of Eq. (6) is 
identical with the right-hand side of Eq. (2); it is shown 
below, for the special case of a rectangular section under 
a two-dimensional temperature distribution, that this 
term is usually much larger than all the others and 
that use of Eq. (2) may therefore be considered per- 


missible in practice. 


DEFLECTIONS OF RECTANGULAR BEAMS 


The shear stresses in a simply supported rectangular 
bar of depth 2c are 


TyrD = [(c? — 4°) ‘2T] (dM p dx) ) 


(S 
Trp = 0 f 
Then 
(Ov/Oyv) = 0; (0¢/O0s) = (c?/2) (dM p/dxo) 
(Sa) 
Substitution into Eq. (6) gives 
vA 2 
v(x) = (1/ET) | | | l(@ET — voyr) yMpn — 
J0 WA 
(1 + v)y?rz,7(d.Wp/dxo) |dxody dz (Sb) 
provided that use is made of the last of Eqs. (4). The 


last term of Eq. (Sb) may be simplified as follows: 


el . 
| f | Vtryt(dM p/dxe) dxody dz = 
0 A 
be OMprryr OTryT 
f | | | ee as Mp = "| dxo dy dz = 
0 Ae E OXo re) Xo 


xo= I Zz, . . 
f j y*( M ptryr) dy dz + / | Mp y° x 
A xo = 0 J0 J Ae 


Ooyr 
dxdy dz (9) 
y 


[CAL 


'?. In what follows, only the deflection v 
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where use has been made of the equilibrium equatio 


(Orryr/OxX) + (Oo,7r/OYV () a 


All surfaces are, however, free of tractions; then 7 


zero when x) = O, Z and, therefore, 


v/ . 
u(x (1/EI) | | 
O e fe 


: Od,7 ‘ : 
(l +p cz Mop dy.dy dz Qh 
Oy | : 


For the present case a, 
integration by parts gives 


g . Oo V7 s 
y dy = (y*ay,r) —2 vy o,rdy (9% 
Ov ‘ 


where the first term of the right-hand side is zero 
Substitution of Eq. (9c) into Eq. (9b) finally gives the 
deflection as 


| o/ . . 
a(x | [eET + (2 + 
FI . 0 le 


V)Or7 |.\Tpyv x 
dxydy dz l() 


Comparison of Eqs. (10) and (2) shows that the 


only term omitted in the latter equation is that pro 
portional to the quantity 
*/ 
| 
0 


Oy7 Mop 7 dxody 


where 6 is the width of the beam in the z direction 


Reference 3 gives, for a statically determinate rec- 


tangular beam subjected to an arbitrary two-dimen 
sional temperature distribution, the following result 


o,7/ak = ‘ —y | T"dy + { T"y dy + 


v 
(c/4) | T” dy (1 + 2(y/c) + (v/c 


e 


(1/4) | T"y dy({2 + 3 (y/c) 


1c 6) | T’’" dy — (y*/2) | TT" y dy 


« 


(y/2) | T’'"’ y? dy — (1/6) , ef 


2(y/c)* + 2(y/c)*] + 


(c¥/24) | T’’”’ dy{(y/c)? + 


(c?/40) | T’’” ydy X 


[4 (y/c) + 10 (y/c)? + 7(y/c)* — (y/c)*] — 


(c yf T’’"’ y? dy [1 + 2(y/c) + (y/c)?] + 


(1/24) [ T’’”’ y® dy [2 + 3(y/c) — (y oy r 


i} terms containing only sixth and higher derivatives 


of T! (il 


must vanish at y +c; then 


we 
thi 
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CALCULATION OF THERMOELASTIC BEAN DEFLECT 
: 
quation where primes indicate differentiation with respect to TABLE 1 
) Inspection of Eq. (11) shows that, for temperature Per Cent Error Introduced by the Use of Eq. (2) in Place of 
d) § ; . : E 10 
f distributions which are either constant along the span ’ 
s 4 } ” ” = ) ” ) 
len 7,7 1s or varv linearly with x, the stress o,7 1s zero; Eq. (2 
‘ . 2 7) Oo” ()&7 (ye 
is therefore exact in these cases. For other temper il Qn — 
| ature distributions Eq. (2) will be in error, but if the 2 0.003 0.5 0.4 
' . a , =a 3 0.000 cs 0.9 
variation of 7 along the span is sufficiently smooth the Ff 0 09 19 15 
higher derivatives of 7 with respect to x will be small, 5 ).6 2.8 2.2 
j oy } the term (2 + v)o,r will be small compared with a7, 
IZ 18) 5 . . ° . 
and therefore Eq. (2) will provide a good approxima calculated separately from Eq. (1) and then superim 
be th tion. This qualitative discussion may be made more posed on those just obtained from Eq. (2).* 
=—C, el - : a - - ; s as . - é 
specific with the aid of the following example. It is of interest to note that the accuracy of Eq. (1) 
Let the deflection in the y-direction due to a tem may also be examined in a manner entirely analogous 
/ perature distribution of the form to that just employed here for Eq. (2), by using the re 
ay Ye ee py” . i : 
‘ t T ’ ; ee 9 sults of reference 4 in place of those of reference 3. It 
1 To(x/L)"(y/e (12 ; . : 3 ; 
is shown in reference 4+ that an expansion of the same 
po be calculated at the center (x L/2) of a simply sup- form as Eq. (11) may be written for the stresses due 
Sives t te ean The > > ic > 4: . . 
ws " ported beam. Then the moment 7p is to any type of distributed loadings on a rectangular 
: , , beam 
fxoll — (x/L) Por Ss Xs 2 oe 
; S > - ' ; 
Mp(Xv, A = ae ee ee , 1s In conclusion, it should be remarked that in the 
— x[(x0/L li forx < % <1 ; ee 
} analysis presented above the effect of inertia was tacitly 
. . . - “i J - 
dz (10 and the deflection may be written in the form neglected. The error thus incurred was examined in 
I o(L/2)/aTob 2L?2 = v9 + (c/L)2m + reference 6, and was found to be small in usual engi 
| o(L/= / u 4) 0 : eee 
that the neering applications. 
hat pro- (c/L)iw 4 (14 
. : ' ; a ; REFERENCES 
Here w% is the contribution of the term aT in Eq. 
10) and, therefore, comprises the entire deflection as ley, Bruno A., and Mechanic, Harold, Therm 
i : i 9 a and Deflections in Beam-Columns, WADC Tech. Rep. 54-425, 
calculated from the approximate Eq. (2). The suc ine 
: . > March, 1954 
—e ceeding terms are the required correction due to the 2 Timoshenko, S.. and Goodier, 1. N., Theory of Elasticits 
-CT1O : 7 gee » »., al x r, J . q y lasticit: 
nix ald o,r term; v; results from the first bracket of Eq. (11), Second Ed., McGraw-Hill Book Co., New York, 1951 
»-din v» from the second, and so forth. It may be seen from ‘Boley, Bruno A., The Determination of Temperature, Stres 
~aimen- . 7 eff s i wo ensiona rermoelastic Problem 
of Eq. (14) that the shallower the beam the better the end Deflecttons in Two-Dimensional I hermociasi Protiem 
result : - : s : Journal of the Aeronautical Sciences, Vol. 23, No. 1, pp. 67-75, 
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- an upper limit to the error to be expected in practice. * Boley, Bruno A., and Moore, Roland H., Note on the Calcu 
ue ; es lation of Deflections of Indeterminate Structures, Journal of the 
The results are collected in Table 1; the maximum ; ci es cae! meal . 
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error encountered in the range of exponents considered See also Vol. 18, No. 4, p. 285, April, 1951 
is only about 2.8 per cent of the correct value. It is 6 Boley, Bruno A., Thermally Induced Vibrations of Beams, 
| therefore reasonable to conclude that Eq. (2) will be Journal of the Aeronautical Sciences, Vol. 23, No. 2, pp. 179 
y adequate for most practical calculations. 181, February, 1956 
The above developments have been limited to stat- * Note that in the calculation of the thermal deflections of 
| ically determinate beams, but the conclusion just statically indeterminate beams the simplification discussed in 
Bg : _ : . : reference 5 can still be employed; in other words, it is permis 
stated is easily extended to indeterminate beams as ; ; ne 
a ine “ se ; 2 . sible to disregard entirely all redundant quantities (and there 
i well. n fact, the additional deflections occurring in fore to treat the beam as a determinate one) in the calculation 
y+ 4 this case, due to the various redundant loads, may be of the stresses due to the dummy loads 
y X 
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Optimum Thickness and Lift Distribution for 
a Fuselage in the Presence of a Prescribed 
Wing System 


M. E. GRAHAM* 


Douglas Aircraft Company, Ine, 


SUMMARY 


\ drag formula of Hayes is utilized to obtain a generalized 
“area rule’’ expressed in terms of sources, lift, and 
Simple examples of optimum fuselage 


supersonic 
side-force singularities 
source distributions and optimum fuselage lift distributions are 


presented. 


DISCUSSION 


§ hes RULE, in supersonic flow, for the optimization 
of a fuselage-thickness distribution in the presence 
of a prescribed wing-thickness distribution has been 
described by various authors—e.g., Ward! and Jones.° 
It is interesting to see how easily this rule is derived in 
terms of singularity distributions from the drag formula 
developed by Hayes,* and how, from Hayes’ formula, 
the generalization of the rule to include lift and side- 
force as well as source distributions is immediately 
evident. Such a generalization, in terms of area dis- 
tributions and forces, has already been presented by 
Lomax and Heaslet.4 However, the subject seems 
of sufficient interest to be reconsidered from a slightly 
different viewpoint. Since the drag formula of Hayes, 
which is valid within linear theory, is expressed in 
terms of singularities, it seems natural to determine 
optimum configurations in terms of singularities and 
then later obtain the geometry of the configuration. 
(The analysis presented here appears also in a slightly 
different form and in somewhat more detail in refer- 
ence 5.) 

For midwing monoplanes, the rule in its usual form 
(without forces) is adequate to prescribe the optimum 
fuselage-thickness distribution. However, for con- 
figurations such as high-wing monoplanes and ring 
wings, the thickness distributions and the lift and side- 
force distributions do not act independently in produc- 
ing drag. Thus, in general, the optimization of the 
fuselage-thickness distribution is affected by the lift 
and side-force distributions on the wing as well as by 
the wing thickness. In addition, there is the possibility 
of optimizing the fuselage lift distribution. 

This generalized rule prescribes what lineal dis 
tribution of sources, lift, and side force should be used 
to represent a fuselage in order to minimize the wave 
drag of the fuselage in combination with a given wing 
(the wing being represented by a spatial distribution 
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of singularities). Furthermore, if no net lift or sid 
force is carried on the fuselage, the minimum pressure 
(No net 


lift or side force on the fuselage means that the fuselag 


drag (wave plus vortex drag) is achieved. 


singularities, which are distributed along a streamwise 
line, will produce zero rather than infinite vortex drag 

Although the actual geometry of the wing-body com 
bination has not been obtained, the fuselage produced 
by these singularities is cambered and has essentially 
circular cross sections. Distortions of the fuselage 
cross sections, which are considered by 
Heaslet,* have not been included in this analysis 


Within slender body theory the rate of change of the 


by Lomax and 


fuselage cross-sectional area is proportional to the 
strength of the fuselage source distribution. 

A wing fuselage combination may be represented by 
disposing singularities along a streamwise line, say the 
é-axis, to form the fuselage, and by distributing singu 
larities in a region of space to form the wing system 


The wave drag of the combination is, according t 


Hayes,* 
l (27 -¢gffO Bsin 6 | 
Dp = | 1 (&) — hig 
2a Jo 62a J OF: 2g 
B cos 6 Des 8B sin 0 
f,(é1) + WE, 6) - (£6 
2q 2q 
B cos 6 W Oo 8 sin @ 
r(£1, 8) i Ags) — g 
Dg F 0% 2g 
8 cos 0 ae B sin 6 
+ — Sale) + Wsl&, 6) — — V")(0, 0 
2q 2q 
Bcosé@__. 7 x ; —? 
IV, (&, 0) | In| &y — & dé déod0 
2q 


where the inner integrations extend over all the singu 
The subscripts s, /, 4 denote sources, lifting 
f(&), where 7 = 5, 


larities. 
elements, and side-force elements. 
/, or h, are the fuselage singularity distributions 
W(é, 0) are “equivalent” lineal distributions formed 
for each radial angle 6, measured from the horizontal, 
by transferring the spatial wing singularities w,, w, 
and w, along Mach planes to the £-axis. 

This equation (discussed in detail in reference 6) was 
derived from linear theory by considering the momen- 
tum flux at a distant control surface enclosing the spa- 
tial distribution of sources, lifting elements, and side- 


force elements. 


It gives the drag of a spatial distribu- 
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tion of singularities in terms of the drag of an equiva- 
lent ‘lineal’ distribution. 

It is important to note that, since the fuselage singu- 
larity distributions f,, f, and f;, are actually lineal, they 
After interchanging the order 


do not depend upon @. 


- 


of the @ and é integrations, it can be seen that the fuse- 
lage singularities interfere only with certain averages 
To illustrate 


combination 


of the wing singularity distributions. 


this, consider a wing-fuselage formed 


only from sources. The wave drag becomes (omitting 


the subscript s 


-@ « / O a oO vita ' 
D ; } f(&1) I\S2) + 
med | (Og, ° O&» 
o oe 
° f (£1) . W (go, 0)d6é + 
Of; Of 2/0 
) Oo 1 
f(t) | WW(&, 8) do' x 
a" Bh Bete 
| vIn ann f . > O 
ines — Gldede + Lf "IAP 2 x 
27. 0 ZT. * Of; 
O 7 , 
W(é;, 6 W(, 0) In & — &) d&id&d@ (2) 
0£: 


Now the fuselage 


ve @ of the wing sources where @ = |] 


ag 


Ws, 0)d@. If @, which is a lineal distribution, is com 
bined with f, then the first two terms can be rearranged 
to form D,,,, the drag of the combined fuselage and 
erage wing distributions, and D;, the drag of the 


werage wing distributions. 


D foe. | 2a SS (o Of; [f(g — W( &; (O O£>) x 
[ f(g) + @(k)] In | & — bo) d&déo 
D g 2a) SS (O OE )D(E) (0 DE) DE) X 
In | & £& d&dé 


lhe final term of Eq. (2) is just the drag D,, of the wing 


alone. The drag of the wing fuselage combination 


thus becomes 


D.., = Ds + D, (3 


lf the wing system is fixed, then only the first term 


of Eq. (3) can be minimized. Thus, the fuselage source 


distribution f = (f + @) — ®is an optimum if (f+ @ 
isan optimum source distribution. 

If lift and side-force elements as well as sources are 
included, the derivation is analogous to that above, and 
Eq. (3) still results. However, in the equations for 
D..,, and Dz, f is replaced by > y¥,(0) f(E) and® 


; B 
by }° y,(0)@,(¢)—where y,(0) = 1, y,(0) = — 
hu Ms 
-q 
sin 6, and y,(@) = —(8/2q) cos 6—and the drag ex- 


Note 


now consists of lift and side-force elements @, 


that Z@ 
and @, 


expressions are each averaged over @. 
as well as sources @,. Each of these @; can depend 
upon all of the different types of wing singularities 
i is 


W., W), W,. (See Example 1.) The definition of @ 


x 7 (0) IV; (é, 0)dé 


o2s 
7 :7(0) dé 
0 


(The optimum strengths of the three lowest ordered 


multipoles of Lomax and Heaslet* correspond to the 
negatives of these @. 

The generalized area rule thus states that the opti 
obtained 


mum fuselage singularity distributions are 


if (f + @) consists of the optimum source distribution 
(f; + @,), the optimum lift distribution (f/f, + @,), and 
the optimum side-force distribution (f, + @,). [(/-4 
@) can be separated in this way since a source and a 
lifting element (or two other unlike singularities) lying 
on the same streamwise line have no interference drag. 


For certain important cases (f + @),,,is already known. 


Example 1—Optimum Fuselage Source Distribution 


Consider the problem of optimizing the fuselage 
source distribution f, in the presence of a given wing dis 
Let the fixed 


, and require the fuselage to have the 


tribution of sources, lift, and side force. 
wing volume be V 
V;. 


@, can be pictured as producing a 


fixed volume The lineal source distribution f 
fictitious body of 


revolution of fixed volume 


| az | f.(E1) + @.(E1) |dé 
. 


(It is sufficient, though perhaps not necessary, to assume 
slender body theory in order to say that V,; + J, is 
the physical volume produced by the wing and fuselage 
sources. This is discussed, in part, in reference 7. It 
is not the purpose of this paper to make distinctions be 
tween the use of slender body theory and linear theory 
to determine the geometrical configuration produced 
by the singularities. ) 

If the wing and fuselage lie in the same horizontal 
plane, it can be shown that the wing lift distribution 
does not contribute to @,(£).. Presumably such a wing 
carries no side-force distribution, so @,(&) depends only 
Then V Vy + V 
the case of the midwing monoplane for which drag due 


on the wing sources. This is 
to thickness is independent of drag due to lift. It is 
to this case that the usual supersonic rule for optimizing 
the fuselage-thickness distribution applies. 

However, in general, @,(£) does depend also upon the 
Then V’ V,- 
dependent 


wing lift and side-force distributions. 
V. + V;, 
the 
monoplane or a ring wing with 


where |’, is a “pseudo volume”’ 


upon forces on the wing. A Iifting high-wing 


an outwardly directed 
a negative 


normal force distribution each contribute 


pseudo volume. (For a discussion of pseudo volumes, 
see, for example, reference S. 

If the length c, of the fuselage is specified, and if the 
fuselage extends at least over the mean lineal wing 
distributions (i.e., the wing system lies within the double 
Mach cone defined by the fuselage), then the fictitious 
Thus, 


body produced by f, + @, has the length cy. 
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under the constraint of fixed length and fixed total REFERENCES 
volume, (f, + @s)o; iS a source distribution obtained 1 Ward, G. N., The Drag of Source Distributions in Lincari- 
by Sears® and Haack.'’ (In a practical case, other Supersonic Flow, The College of Acronautics, Cranfield, Report } 
constraints, such as minimum admissible cross-sectional No. 88, February, 1955 
area for the fuselage, might be necessary.) 2 Jones, R. T., Theory of Wing-Body Drag at Supersoy 

Speeds, NACA RM A53H 18a, September, 1953 

Example 2—Optimum Fuselage Lift Distribution Hayes, Wallace D., Linearized Supersonic Flow, Nort! 


: ke . ” Rgterinan Aviat a a See ane 
If the net lift, L,, carried on the fuselage and L,, car- maine TS, SHE. Lae Sapnes, Rapes. ie, AL. 
J June, 1947 


ried on the wing are fixed, then the lineal lift dis- 


* Lomax, Harvard, and Heaslet, Max. A., Recent Developmen 


. 
tribution f, + @, produces a fixed lift L = | in the Theory of Wing-Body Wave Drag, presented at 24th Annual 
Meeting of the Institute of the Aeronautical Sciences, New York 


(fi + @)d§ = Ly + Ly. If the fuselage has a specified January 23-26, 1956. (Preprint No. 617 

length c, extending at least over the mean lineal wing ’ Graham, M. E., Application of Drag Reduction Method 
distributions, then the lineal distribution f, + @, also Supersonic Biplanes, Douglas Aircraft Company Report 

has the length c,. Thus, under the constraint of fixed SM-19258, September, 1955 } 


length and fixed total lift, the optimum distribution for *Graham, E. W., Lagerstrom, P. A., Licher, R. M., and 
Beane, B. J., A Theoretical Investigation of the Drag of Generalized 


ies teas aie ‘ Be an Aircraft Configurations in Supersonic Flow, Douglas Aircraft 
it is likely that additional constraints, such as limiting Company Report No. SM-19181, July, 1953 . 
- 7) an) epo NO. Si PII, | , Yoo 


the magnitude of the local angles of attack of the 


f(&) + @,(&) is elliptic.'' (Again, in a practical case 


’ Bleviss, Z. O., Some Integrated Volume Properties in Line 





fuselage, would be comreniaiabes ; ised Flow and Their Connection with Drag Reduction at Supe 
If no net lift is carried on the fuselage (L, = 0), then sonic Speeds, Douglas Aircraft Company Report No. SM-19469, [ 
the fuselage produces no vortex drag. Eq. (3) then February, 1956 
applies to wave plus vortex drag, if the term D, is ‘Graham, E. W., Examples of Wave Drag Calculation by At vi 
interpreted as the wave plus the vortex drag of the wing. plication of Fourier Analysis about a Streamwise Axis, Douglas 
It is interesting to note that Hayes’ equation pro- Aircraft Company Report No. SM-19529, January 1956 A 
vides a certain equivalence between lift distributions * Sears, W. R., On Projectiles of Minimum Wave Drag, Quar- | 
e ° e a ° ° ° ° wr rT se. ‘ mM: CS ] r ‘ a 9 7 
and source distributions. The lift distribution f)(£) + terly of Applied Mathematics, Vol. 1V, No. 4, January, 194 
- : ~ ; : . 0 *K y M osstorinie qi Sle ‘ellenwi erst l¢ ' 
@,(&) can be converted into a source distribution which " Haack, W., Geschossformen Kleinsten Wellenwiderstand Ap| 
= : ae — Bericht 139 der Lilienthal-Gesellschaft fiir Luftfahrt, 1941 
would produce the same wave drag. This source dis- = aint awa as 
‘buti i cae : : \lso available as Translation CGD-253, Goodyear Aircraft 
tribution represents a fictitious body of revolution with Corp., or Translation No. A9-T-3, Brown University, 1948 * 
a base area prop al to the tot . The imu ae : 
i b se area Proj ortional to the t at ullift. The optimum 1 Jones, R. T., Theoretical Determination of the Minimum Drag 
elliptic lift distribution mentioned above converts of Airfoils at Supersonic Speeds, Journal of the Aeronautical a 
(within slender body theory) into a Karman ogive. Sciences, Vol. 19, No. 12, December, 1952 : 
s : g : 
} oe 
Eff 
s 
rr . . . . ’ . . y* Tot 
A Two-Dimensional Approximation to the Unsteady Aerodynamics of Rotary Wings 
Continued from page 92) set 
4 
REFERENCES § Schwarz, L., Berechnung derdruchverteilung siner harmonisch Th 
H j he D sich verformenden tragfiache in ebener stroémung Luftfahrtforschung, ( 
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Annual Forum, American Helicopter Society, Washington, D.C., © Bierens De Haan, Nouvelles tables d’integrales definies, Stech 
May 14-17, 1953. ert and Co., 1939. 
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* Isaacs, R., Airfoil Theory for Flows of Variable Velocity, 1! Séhngen, H., Wathematische seitschrift, Band 45, pp. 245 I 
Journal of the Aeronautical Sciences, Vol. 12, No. 1, pp. 113-117, 264, 1939. 
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_ om . , . < > > 7 r 4 , of Ae na nile stadiuity 
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Vibration of Clamped Parallelogrammic 
Isotropic Flat Plates 


M. Hasegawa 
Shimane University, Matsue 
September 6, 1956 


Japar 


V IBRATION of rectangular plates was discussed by W. Ritz,! 
S. Iguchi,? D. Young,’ and M. \ 


skew cantilever plates was discussed by M. \ 


Barton,* and vibration of 
Barton Vibra 
tion of clamped rhombic plate was discussed by the author In 
of the 


tion of clamped parallelogrammic isotropic plates by the 


this paper, we shall calculate the smallest frequency vibra 

Ritz 
method 

Let the oblique coordinates axes (x, y) be taken in the middle 
plane of the plate of uniform thin thickness in such a way that 
the origin coincides with the center of the plate and the axes art 
parallel to the sides. Let us denote the lengths of side, thickness, 
the density, Young’s modulus, and Poisson's ratio of the plate 


P, E 


coordinates axes be denoted by a If x 


by 2a, 2b, h, and », respectively. Let the angle between 


is the transverse dis 


placement of a point on the middle plane and ¢ the time, the dif 
ferential equation for w is 
(D/sin* a) [(0*w/dx*) — 4 cos a(dtw/Ox dy) 4 
2(1 + 2 cos? a) (0*w/Ox*0y"?) — 4 cos a (O*w/OxXO4 
ph(d*%w/ Ot?) = O l 

where D Eh? /12(1 — v? 
The boundary conditions are 
w= 0 and Ow/dn = 0 at x = +a and ¥y th (2 
where 0/0n denotes differentiation along the normal 

\s the problem is so complicated, we shall use the Ritz 


method According to this method, the mode of vibration is as 


sumed If we assume the following mode, 


w = W(x, y) cos of 

Eq. (1) is transformed into 
(O4NW'/Ox4) — 4m(O*W/dx80y) + 201 + 2m?) (O4W/Ox*0y?) — 

$m(O4W'/Oxdy*) + (O1W/oyv') — A’W=0 (4 
with the boundary conditions 
W=0 and OW/on 0) at 2 ta and y = +i 
in which m = cos a and 

K’' = (ph p* sin* a)/D 6 


In order to apply the Ritz method, we consider th 


integral: 





D i *b orn 

2U = - dx (AH)? — 4 cos All 
sin? @ | , F . 2. Oxoy 

4 cos? a(O071/dx Oy)? + 2(1 — v) sin? a) (0711 /Oxdy)? 


(221 /dx?)(02W/dy?)} — K | dy (7) 
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We assume for the expression of 1’ which satisfies the boundary Ol'/Oa 0; OU/da 0; Ol'/da 0) 0L'/oa } A 
conditions of Eq. (5) that OU/0ax2 = 0; Ol’/da 0; OU'/da 
W = (a? — x?)? (6b? — y?)? (aon + auxy + doox? + aooy? 4+ hese equations are linear in the A’s and are found to be j 
l pa y3 22 8 
ee ee ee BywoAw + BurAu + BioA%® + Bide + Bin A 
in Which doo, @i1, and ds are arbitrary constants and Eq. (8 BysAi3 + B I U 
is introduced into / 
Since the constants are to be determined so as to make L’ a By Ao + ByyAn + BrwA + BrorAg + BryA 
minimum by the Ritz method, we must have Bry3A + BrwA 0 
4 
- — - _ where 4 
ee } 
b =A pha‘p? sin’ a 
| = bi*iq.., k= i K'a‘ = 
‘ a D 
: l+k 1k2(1 + 2m? 2Kk 2m(k? + k 
Bio ™ — +4 _ B a 
1,575 11,025 Y9 225 11,025 
k l 2Kk* 
a, \ ee 
11,025 17,325 1,091,475 
7 k k ki + ] $k2(1 + 2m?) 2Kk 
B = B 2m t . ££ : 
33,078 10,425 24,255 99,225 12,006,225 
: ’ : ki + 1 2Kk4 
B = 8 = ¢, = hb; = _ B = 
121,275 12,006,225 
: k k 2m(k* +k 
B = Bs, = 2m oe — i By, = = 
33,075 121,275 363,825 
k l 8k2(1 + 2m? 2Kki 
B = B = r | . = B, 
72,465 105,105 1,091,475 52,026,975 
ki ] $k2(1 + 2m? 2Kk 
b = —+- | —_— = 8 
3,679 79,075 121,275 1,729,725 j 
; ks k k k Fu 
Bb Bs = 2m } , B = Bio = 2m — =f 
121,275 525,525 121,275 363,825 
: A 2Kk ‘ 
5 = B B = B 
1), 425 525,525 52,026,975 
13k l 68k2(1 + 2m?) 2Khk 
B = { | : ; = B 
800,415 315,315 14,189,175 156,080,925 
; kt +1 16k2(1 + 2m?) 2Kk ; | 
B Bex, = — — of -, B ee 
315,315 12,006,225 225,450,225 i 
| k k | 
I = By, = 2m {- ——< ee = B; 
1,334,025 1.5/6,575 
kt 13 68k2(1 + 2m? 2Kk* 
B, = 4 1 ome 
315,315 800,415 14,189,175 156,080,925 
ki + 1 $k2(1 + 2m’) 2Kk* 
B } 
175,175 1,334,025 225,450,225 fi 
The coefficients By, By3, Bao, Bos, Bas, Bas, Bye, Bay, and determinantal equation are calculated for several values of « 
Bee are obtained, respectively, by interchanging the numerators k and are shown in Table 1 
of the first two terms of each of the following coefficients: By»29, The exact value of A for the square plate is 80.9348 which i 
Bix, Booo, Box, B32, B31, Bs22, B33, and B For instance, we obtained by Prof. S. Tomotika.* We shall report the more ac- 
have curate values in a subsequent paper | 
Byo = (1/11,025) + (24/17,325) — (2KR1/1,091,475 —— 
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“ " A Useful Method for Integrating the This curve varies with x between 1 and 0 if we set A B 





9 Time-Dependent, Viscous Equations of 1/2. As ¢ increases, it is displaced parallel to itself without 
to be Motion* deformation. The slope, alwavs negative, reaches the extremut 
} - is 
B?/y, which can be made arbitrarily large for proper 
' Filler. H. F. Ludloff, and P. Treuenfels Therefore, it may be considered a propagating shock in a vi 
) l 
New k University, New York cous medium 
: 195¢ Using an implicit scheme, the difference equation i 
() u;, u;, 1 os 
s A PRELIMINARY to numerical integration of the three : > + 4 
. : : ; ; At 2. Ax 
t dependent, viscous equations of motion, it is in 
: WMAE : ) 
structi to study the ‘‘Burgers equation,’’ which contains bot] _ ’ att 
iracteristic features of the motion of a real fluid—viz., the An 
well as the viscosity effects rhe most interesting 
if ce procedure to be used for integrating this parabolic Here indicates the space point and & the time cycl Using 
s . . - n’ ] voc ’ } } 
differential equation (as well as the true equations of motion von Neumann’s error analysis, it is easily seen that the 
@ ° , ; , 1c c cf} fn / oO ; and 4 
is an implicit scheme, in which the time steps may be taken condition is satisfied : 1 Af) Av and v/ An 
x y large We have carried out an analvsis of stability In order to solve the implicit difference equation | Eq 
id of convergence of iterations before carrying out the nu iteration procedure, indicated below, has been employed. It 
ierical integration and comparing it with an analytic solution ippears appropriate to use here also von Neumann Fourier 
ich we had derived before inalysis to derive a criterion for the cor ence of the iteratio 
An analvtie solution of procedure We use superscripts to designate the number of 
iterations carried out in a time cvck It is best to collect t 
it > e = l it . " 
two linear terms at lattice point 7, k + 1, in iteration cyck 
t 1+ B-tanh[(AB-t B-x)/21 n + 1, while all other terms are kept in iteration cycle ». Ther 
t ¢-(At/2An t u v/ AX Mt An d ( 
1 + (2»/Ax At, Ax 
We substitute u = ij, it 6 
} where 1;, is the exact solution of Eq. (1) and 6 , the iteration error committed up to the wth iteration of the (k+1)th time cyck 
Furthermore, note that “;, , = a@,, , without error since iterations in &th cycle were carried to completion Phen 
j I 
; u-(At/2Ax) (6 6 1) t+ (v/Ax At, Ax 6 + § 
' 0 
| 1 + (2v/Ax At/ Ax 
After introducing 6), .+, = e'7?>*-e?", we obtain for convergence 
u-(At/Ax) (1 sin BAx) + (2v/Ax At/Ax)-cos BAx 
; 2y/Ax)+(At/Aa 
a es _ — — 
= C,? + (C.? C,2)-sin? (8-Ax) < 1 - 
an At/Ax Sv 1/[u 2v/Ax 
wher — . . a . 
This convergence criterion has been checked for the values 
2v/ Ax At/ Ax ; u-(At/ Ax of At/Ax and v/Ax, indicated by crosses in Fig. 1 (where A 
rf : r( so a 7 
1 + (2y/Ax)-(At/Ax 1 + (2y/Ax At/ Ax B = 1/2), and checked perfectly with difference calculations 
p , — . ee F : carried out in the N.Y.U. thesis of L. Goldberg, divergence of 
If S. 2y/Ax, then C2 \ GC, and | e’| mar? = Ci? NS 1 is satis , ; ee : ; 
. : nig : — : iterations resulting in increasing oscillations Also, the nu 
fied for any values of At/Ax If «w > 2v/Ax, then ( > Cj, and : j ; sd ; 
BF eed al ‘ merical shock curve advanced in agreement with the analytica 
C2? S 1 vields 9 
ia curve, although for approx. At/Ax > 2, gradually changes of 
shape occur—an indication that the solution of the differenc« 
Investigation sponsored by ONR and ORS equation gradually deviates from that of the differential equa 
OF @ an tion with increasing mesh width 
XY Analogous considerations have been carried out for the sys 
whe 3 . ° . - 
which is tem of the three viscous equations of motion For this pur 


more ac- pose, new concepts have to be introduced which will be explained 


elsewhere 

















7 Approximate Trajectories of High-Drag Bodies 

nive t\ 

J. Api George S. Campbell 

NJ Research Engineer, Hughes Aircraft Company, Culver City, Calif 

> 1 
, September 10, 1956 

Pi 

a € 8 O : 4 INTRODUCTION 


kyo Im | u AX § ions NOTE presents charts for the rapid determination of 
Fic. l. approximate trajectories of high-drag bodies. Such a prob 
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50¢ initial deceleration, 49 = O 


lem might arise, for example, in the retrieval of a missile by 
parachute deployment 
illustrated in Fig. 1, in which Dp represents the drag of a para- 
chute or similar device. Within the accuracy of this analysis, 
the actual source of the drag is immaterial, and so only the total 
drag D is considered. Furthermore, the lift acting on the 


niissile is neglected 


(NALYSIS 
The tangential and normal force equations are 
mV = —D + W sin 4 la 
mVy = W cos ¥ lb) 
Although these equations are relatively simple, it is necessary 
to make further assumptions before a closed-form solution may 
be obtained \ useful simplification results if W sin y is re 
placed by HW’, which is the asymptotic value of the weight term 
The physical justification of this approximation depends on the 
assumption of a high-drag body—i.e., D > W, initially. The 
approximate solution obtained in this manner may be regarded 
as the first step of an iterative procedure 
First Approximation 
The differential equations for the first approximation to the 
solution of Eq. (1) are 
Vi + AV? = g 2a 
Vivi = g cos 7; 2b) 


rime histories of velocity and flight-path angle are obtained by 


simple integration of Eq. (2) so that 


Vi/Vo = Betnh u | 
= 2tan~' (Ct) — (4 2) f 
The variables u and é are related to the time ¢ by the equations 
u= A+ BRVot 
— = cosh u 


and the nondimensional constants 4, B, and C depend on the 
initial deceleration and flight-path angle: 


A =tanh'B 


The forces acting on such a vehicle are 


C = tan [(y0/2) + (2/4)]/cosh A 


Altitude loss and horizontal distance from an initial p 
are obtained by integrating the proper components of the 


ant velocity ; therefore, 


Kith ho) = 


tan ct tan Cz Cetnh™! & + Cctnh 
Second Approximation 


In many cases, the trajectories calculated from Eqs 


+) may be of sufficient accuracy for engineering use. How 
if more accurate results are required, the following corr 
procedure may be used The main source of error in th 


approximation is the use of a constant weight term in the 
equation, which results in predicted velocities that are too large 
rhe difference, Al’, between the true velocity and the fir 
proximation may be obtained from the linear differential 


tion 
d/dt) (AI + ZA V (AV) = g(sin > 


in which J and are given by Eq. (3). This relation is derive 
by subtracting Eq. (2a) from Eq. (la) and linearizing the result 


ant expression. The solution of Eq. (5) is 


AV B VC ly 
sinh? x Cc? a 


C? + cosh? u +2 VY C2 + Lleosh uw sinh aw + (C? + 1) sinh 
C? + cosh? A +2VY C?2+ 1 cosh A sinh A + (C? + 1) sinh? A 
Vv ¢ l l C? cos! 2 cosh sinh 1 
log log 
( i * cos ! ( cosh 4 + sinh 1 | 


where the approximate value for B has been used 


o> ~ ry 


Bb? ~g KI 





rhe increments in the horizontal and vertical positions du 





the velocity correction, AV, are readily calculated by graphi 


integration 

















sh= f(a sin dt} 
Av = i AV) cos dt\ 
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(using Eq. 4 
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2.8 } Let us begin with the simple, divergent channe r diffuser 


< 


RADIANS one of the basic systems of fluid mechanics It is well kn 





that diffusers are inefficient when the divergence angle is 


oOo00 Oo 
—“NMNWLDO - 





LO 














diffuser with large divergence angle, the efficiency can actuall 





great In fact, it can be shown that, for a two-dimensional 


ta be improved by insertion of a plate along the center streamline 
O — 1% see Fig. 1] This improvement in efficiency results even thou 
= friction surface has been added. The usual explanation 
a A - the influence of divergence angle on efficiency is well put in th 
> a following quotation “As the angle of flare increases beyor 
the optimum value, however, the adverse pressure gradient ass 
r - a « -< aa oe oe a ian ee oes ciated with the deceleration of the fluid becomes big e1 
is p F ee Sa cause a rapid thickening of the boundary layer and separati 
w : from the wall.’’* 
S ; Let us assume that the Reynolds Number is high that 
. 3 boundary layers are consequently thin. Under these 
= = the reduction of flare by the insertion of a plate does not in 
_ t way affect the pressure distribution ig the orig 
rx walls. The adverse pressure gradients are unchanged. S 
Ox ther than the change of pressure distributi t 
a 0 2 4 6 8 O 12 14 16 18 20 22 24 prevented separation 


iOR 





NON-DIMENSIONAL TIME, K Vt Purning to another phenomenon, 


ration and reattachment near the f a thin-plate airfoil 


7 


Approximate time histories for 100g initial deceleration 


using Eq. 4 small angles of attack. Free-streamline solutions of classi 
wake theory do not account for this reattachment In discussit 
this, Prof. A. H. Shapiro suggested in conversation (1952) tl 
RESULTS : . ; : : : 
the phenomenon is probably related to that of jet mixing 
Fig. 2, comparison is made between the approximate tra 


* Hunsaker and Rightmire, / 


tories obtained from Eqs. (4) and (7) and an exact solution : , 
349; McGraw-Hill Book Company, In New York 


f Eq. (1) obtained by the Kutta-Runge method. The body has 
initial deceleration of 50g and zero flight-path inclination 

Use of the simplest approximation, Eq. (4), is justified when a 

osition error of the order of 10 per cent can be tolerated It is ail 

xpected that the error introduced by using Eq. (4) might be 

msiderably smaller than 10 per cent for larger values of the 

itial deceleration or nonzero (but positive) flight-path incli 


ipproximate trajectory into much better agreement with 


nations. Addition of the velocity correction term, Eq. (7), brings Oy, 
exact, step-by-step result " CS) 

In order to facilitate rapid estimation of the trajectory of a \/- a en Se ee — 

igh-drag body, time histories are presented in Figs. 3 and 4 for y 


) and 100g initial decelerations. The approximate equations "> 
Eqs. (4)] were used for these calculations, and so the initial —P, ATE 


leration is related to the initial drag as follows 
g =] A V?7/¢ 


Kg = pCpS/2W 


Variations in density or drag coefficient may, if necessary, be 





iken into account by allowing A to have a stepwise variation ‘ — a : 


ith time 


Discussion of the Influence of Confining 
Surfaces Upon the Boundary Layer in a 
Channel 


sistant Professor, University of Missour 
September 10, 1956 


HE BOUNDARY LAYER in a narrow channel does not behave 
exactly like that on, say, an isolated airfoil. Interactions 
between the boundary laver, the free stream, and the confining 
walls of a channel are now well understood in the case of super 





sonic flow. However, the subsonic case, especially where ap- 


plicable to the passages of turbomachines, seems to lie in a blind 


Bounnary oF 
several other workers in the field of fluid mechanics and have DEAD ~ WATER Reson 


een deemed by them to be of sufficient interest for presentation 


spot. The ideas presented below have been discussed with 


in this Forum Fic. 2 


let us take not f the sey 
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is well known that, for two-dimensional, turbulent mixing be 

tween a jet and a dead region, the mixing zone grows linearly with 
i divergence angle of about 15° (this should be taken only asa rule 
of thumb Now, if we calculate the free-streamline pattern for 
a plate airfoil with separation, we find the dead-water region to 
be roughly wedge shaped with its edge at the point of separation 

If the wedge angle is greater than the angle of turbulent mixing, 
separation will be unaffected 


will knit the free-stream flow back to the surface 


If it is less, presumably the mixing 
\n analytical 
solution, estimating the losses on an airfoil with frictionless sur 
face but with this kind of separation and reattachment, would be 
most weleome Phis would, undoubtedly, be difficult to obtain 
However, in the case of the channel, even the crudest of investi 
gations vields much useful information 

I re 
Let us say that 


Taking account of the above-mentioned phenomenon « 
ittachment, let us return again to the diffuser 
separation has taken place and that the dead-water region satis 
fies the condition of constant pressure. Some distance past the 
separation point the streamlines will have become parailel The 
ingle bounding the dead region will be the angle of divergence of 
the diffuser (see Fig. 2 

If this angle is greater than the angle of mixing, reattachment 
will not take place If it is less, separation will not take place 
Whether or not this last conclusion is true, the angle of sepa 
ration undoubtedly influences the behavior of the boundary layer 
and, where separation has not taken place, the v/rtual angle 
separation offers to be a factor of considerable influence By 
virtual angle of separation is meant the angle at which the sep 
irated flow would leave the surface 7/f separation took place at a 
given point. Any explanation of the behavior of the boundary 
laver must be in terms of factors within the boundary layer or in 
the immediately adjacent portions of the free stream. Confining 
surfaces some distance away cannot in themselves influence the 
boundary laver. They do, however, determine the virtual angle 
of separation which is a factor within the boundary laver itself 
This angle is a point function, variable along the surface. Its 
influence should be of importance even after actual separation 
has occurred, especially when the dimensions of the separated 
region are small 

From these observations it appears that, in any study of the 
boundary layer in a channel, one should take cognizance of this 
distribution and aspect 


angle as well as the surface velocity 


ratio (this last factor influences secondary flow Conversely, 
if two surfaces have the same pattern of velocity variation and 
of virtual-separation-angle variation, one would expect similar 
boundary-layer behavior, irrespective of the actual geometry 
This last statement can be taken as a specification of the condi 
tions for quasi-similitude. Once it becomes established that this 
angle does truly account for the influence of confining walls, it 
possible to generalize experimental measurements 


should be t 
For example 


from a few typical channels to a variety of others 
it should be possible to predict, confidently, the behavior of a 
mixed-flow impeller from studies of a stationary channel modeled 


upon the factors mentioned above 


Two reservations should enter into the evaluation of the 
ideas presented above. First, the validity of virtual separation 
in explaining boundary-laver behavior need not necessarily de 
pend upon the mechanisms proposed here. Secondly, when one 
attempts to compute the angle of separation, he may find, 
say, that at the separation point it is initially zero, gradually 
growing toa finite value. In an effort to bypass such local effects, 
the author has had some success in restricting the analysis to 
overall mass movements of the fluid rather than to local condi 
tions. For example, in dealing with the velocity recovery on the 
surface of a vane in a turning cascade, the width of the separated 
zone in the plane of trailing edges and the distance of this line 
from the separation point were used to compute the angle 
This procedure reduces the virtual-separation concept to nothing 
more than a vardstick for comparison purposes. As such, how 
ever, it still promises to be useful in the study of narrow turbo 


machine passages 


ICAL 
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One cannot but feel that there is some more exact concept 
lving behind that which has been discussed here Now that the 


need has been made apparent it may soon be forthcoming 


Location of Mach Discs and Diamonds in 
Supersonic Air Jets? 


Donald E. Wilcox,* Alexander Weir, Jr J Niche 3 
Roger Dunlap **** 
Aircraft Propulsion Laboratory, University of Michigar 
Ann Arbor, Mich 
September 17, 1956 
SYMBOLS 
1 throat area of nozzle 
i exit area of nozzl 
dD diameter of Mach disc 
D exit diameter of noz 
Va exit Mach Number computed from I 
FP itmospheric pressure 
Pr exhaust pressure 
P upstream stagnation pressure 
‘ distance from nozzle exit to oblique shock- wave intersect 
or distance from nozzle exit to Mach d i n 
ratio of heat capacitie issumed to be i 


INTRODUCTION 


She COMPLICATED shock and rarefaction patterns existing j 
supersonic jets have long been of fundamental interest i 
the study of fluid dynamics. It is well known that an over 
expanded nozzle (exit pressure less than receiver pressure) leads 
to oblique shock waves at the exit of the nozzle which intersect 
and give the familiar ‘‘diamond”’ pattern. If the degree of over 
expansion is great enough, this pattern is modified so as to ter 
minate the oblique shocks with a normal shock—or the s 
called Mach disc configuration. On the other hand, a sufficientl 
underexpanded nozzle will lead to a similar Mach dise because of 
the focusing of compression waves from the jet boundary. At 
extreme degrees of underexpansion, the oblique shocks and their 
markedly 
configuration is appropriately termed a shock bottl 

effort to experi 


reflection from the Mach disc are curved, and the 


Numerous investigators have devoted 


mental and analytical studies of these phenomena. However 
most of these investigations have been restricted to sonic nozzles 
low-pressure ratios (stagnation pressure to r 


Kelber and Jarvis’ have utilized the method 


and relatively 
ceiver pressure) 


of characteristics to analyze the flow from a nozzle operating 


t This research was sponsored by the Office of Ordnance Research, U.S 
Army Discussions of this experimental program with Profs. T. C. Adam 
son, R. B. Morrison, and H. E. Bailey of the Aeronautical Engineering 
Department are gratefully acknowledged 

* Research Assistant, Engineering Research Institute 

** Lecturer in Chemical Engineering and Associate Research Engineer 
Engineering Research Institute 

Instructor in Aeronautical Engineering 
*** Assistant in Research, Engineering Research Institute 
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Fic. 1. Schlieren photograph of supersonic air jet 
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err ' mely high pressures. Such a technique fails to indi humidity owing to the high storage pressure rhe air leaving 
w that t the presence of a Mach dise so that the phenomenological the storage system passed through two dome-control val 
ing siderations for locating the dise are in question Weir® 7 has parallel, through the test nozzle, and then discharged in | 
} reported information on the variation of disc location with pres itmosphere rhe pressure upstream of the test nozzl 
ure ratio for airflow through two-dimensional slits maintained constant by the dome-control valves ‘ 
The influential parameters in determining the Mach disc loca operated by adjusting nitrogen pressure above the diaphragm 
jon il upersonic jet would appear to be the ratio of driving \ny pressure desired from 0 to 600 psig could be maintained 
in pressure to receiver pressure, the nozzle exit Mach Number, and bv this svstem 
the nozzle divergence angle. Information of this type is meager In order to investigate quantitatively the influence of 
nd limited in the range of variation \ccordingly, the initial geometry, four axially svmmetric nozzles with rounded inlet 
and } portion of this research contract was devoted to a systematic throat diameters of 0.250 in., and conical diverging angles of 1 
xperimental investigation of the effect of these variables hese 15°, 20°, and 30 total ang vere designed to b osha 
sults have been used as an aid to the theoretical study which the l-in. pipe downstre was of the control val\ he Mach 
still in progress Phis paper is devoted solely to the experi Number at the exit was 3.5 in all cases, as comput 
ental results isentropic flow relation 
EXPERIMENTAL APPARATUS AND PROCEDURI 1 [ L 1) 2 7 
| \ir was obtained from a 2,000-psia storage system with a total | 7 {1-4 1) /2) VW | 
volume of 170 cu.ft The air used thus had a very low absolute 
Schheren photographs were obtained of the flow issuing from 
these nozzles at upstream stagnation pressures from 55 to 615 
psia in 50-psi increments \ conventional schlieren system with 
a 6-in. field of view was used, the pictures being recorded on 8 
ctic ' by 10-in. film 
In order to investigate the influence of the exit Mach Number 
on the flow, the above-mentioned nozzles were modified \fter 
each series of runs at one Mach Number, the nozzle was short 
xisting in ened until the exit diameter agreed with that given by the are 
iterest j ratio required for the next Mach Number Phe Mach Number 
an over tested were 3.5, 3.0, 2.5, and 2.0 
| 
ire) leads 
intersect 
of over 
is to ter 
the s 
fficientl 
ecause ol 
ary At 
ind their 
and the 
| 
) experi 
lowever 
c nozzles 
© to 
method 
perating 
arch, US. | 
C. Adam ' 
igineering | 
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i 
Engineer i 
| 
P,/P,=3.24 
j ‘ Fic. 2. Effect of pressure ratio on flow pattern. Exit Mach Fic. 3. Effect of exit Mach Number on flow pattern. Pressurs 


= 41.7. Divergence angk 10 


Number = 2.0. Divergence angle = 10 ratio 
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stagnation pressure—10° nozzle 


RESULTS 


A schlieren photograph of the flow issuing from a supersonic 
nozzle is presented in Fig. 1. The jet boundaries as indicated 
by the optical technique and the curved shock waves extending 
from the nozzle exit to the normal shock wave are clearly visible 
Less visible are the weak oblique shock waves in the region be 
tween the nozzle exit and the Mach dise which are due to imper 
fections in the nozzle. 

In Fig. 2, a series of schlieren photographs showing the in 
fluence of the pressure ratio on the flow pattern in the axially 
symmetric jet is presented. Fig. 3 presents four photographs 
of the flow issuing from nozzles with the same divergence angle 
(10°) and pressure ratio (P)/P, = 42) but at exit Mach Numbers 
of 2, 2.5, 3,.and 3.5. 

Measurements of the distance of the Mach disc or shock-wave 
intersection point from the nozzle exit (x) were made from the 
Fig. 4 presents the results of these meas 
The plots for all 
The ratio (x/D,) is plotted versus 


schlieren negatives 
urements in graphic form for the 10° nozzle. 
the angles are very similar. 
pressure ratio (Po/P,) for the 10° divergent angle nozzle in Fig 
4, with exit Mach Number as a parameter. 


DISCUSSION 


As mentioned, the exit Mach Numbers presented were com 
puted assuming isentropic flow. It is realized that there is some 
departure from isentropy in the nozzle, as evidenced by the 
schlieren photographs. 

Figure 4 indicates that the exit Mach Number, in addition to 
the pressure ratio, is a most influential variable in determining the 
Mach diamond or disc location. The effect of the divergence 
angle of the nozzle is relatively small over most of the pressure 


range. 
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Effect of Thermai Resistance of Joints Upon 
Thermal Stresses 


B. E. Gatewood 

Research Coordinator and Research Professor, Air Force Institute 
of Technology, Wright-Patterson AFB, Ohio 

September 24,1956 


perp Griffith and Miltonberger! have considered th 
effect of the thermal resistance of a skin-stringer joint up: 
the thermal stresses produced by a constant heat source on th 
skin |Fig. 1(b)]. By dividing the cross section of the joint [Fig 
I(b)] into 15 sections and using an electronic differential an 
alyzer, they calculated the thermal stresses for several values of 
the parameters for heat transfer to the skin, heat transfer across 
the joint, and skin geometry. Also, Schuh? has solved nu 
merically the equivalent problem [Fig. 1(a)] for no joint therma 
resistance and has constructed curves for the maximum stresses 
in the web against the skin heat-transfer parameter g = 
Lth/kh,)"*, where h is skin heat-transfer coefficient (B.t.u 
hour ft.? F.), & is coefficient of thermal conductivity (B.t.u./hour 
ft. F.), and L and h, are defined in Fig. 1(a). 

The results given in references 1 and 2 indicate that it may be 
possible to obtain an approximate analytical solution for the 
thermal stresses in the multiweb beam or skin-stringer element 
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Fic. 1. Structural element for multiweb beam. 
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case in Which ¢ is approximately 5.0 at m; = 0 and 8.0 at » 
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Fig. 1) for any value of the parameter g and for any value of 
, int thermal resistance If the skin is assumed to have 
iniform temperature from a constant heat source (this neglects 
effect of the web on the skin but is approximately correct 
t kin? heavier than the web), then the skin temperature 


is time, p. is skin density, (¢ is specific heat, and 7. is 
jiabatic skin temperature 
Carslaw and Jaeger,’ the heat-conduction equation for 


stringer or web is 
OT /Ot = (R/ pe 0?7'/ Ox? 2 
If there is no heat loss at x 0 and if /; is the thermal conduct 
icross the joint from the skin to the stringer, then the 
itial and boundary conditions on the stringer or web are 
] 0 at ¢=0, oT/dox = 0 at x = 0 ' . 


OT /Ox h;/kv) [To(1 — e~74') — T] at x L{ 


Phe solution of Eq. (2) with the conditions of Eq. (3) is 


] 1 — C, exp (—g?Wt/t,) cos (qx/L) + 
Soa, exp (—p,2Wt/t,) cos (prx/I 
La ] 
) 0 
Hpi ¢ L2/k pA / HW Ril; / pole L*. . 
m; = h;L/k 
( ” m; COs q — q sin gd), Pp = mM, cot p 
2g? sin p g? — fr7) (px + sin p, cos p 


there is no resistance at the joint or if the skin and web are 


| structure, then m, ind the expressions simplify to 


rT 
gra 


For the case of no bowing of the structure (symmetrical tem- 
rature distribution or external bending restraint) the thermal 


ress in the web is 
ag | 
—aET(x) + (1/A | aET(x)b(x) dx (6 


v) is the width of the web or stringer and A is the total 

i of the skin and web. Fig. 2 shows the maximum thermal 

stresses in the web (x = 0) against the parameter g for several 

values of the joint thermal conductance parameter m;. The 

es are drawn for a skin and web combination of the same 
1aterial with @ A,/A 5.0, where A, is the skin area and A 

is the web area 
It should be noted that Fig. 2 can be drawn for different values 
f the parameter @. Schuh’s curves? for m; = show some 


variation with @ and also a small variation with the ratio h,/h, 


l issumption used to obtain Eq. (1) loses this latter ratio 
h However, for 2h,/h, varying from 1 to 7.5, the curve 
for » in Fig. 2 varies at most by 15 per cent from Schuh’s 


s, With practically no variation for g > 5. Fig. 2 also shows 
that the thermal resistance has little effect on the maximum 
could be used 


Stresses for ¢ 


5 so the Schuh’s results for m; = 
rectly 
rhe caleulations in reference 1 were made for three values of 


parameter g (approximately 1.55, 5.0, 15.5), several values 


including m; = O and o, three values of ¢ in the range of 

to 12, and h,/h, = 2.0. Note that the parameter @ for the 
ise of reference 2 [see Fig. 1(b)] is different at m; = O and m, 

flat flange of the web should be included in the web area 

for QO but not for m,; If the proper ¢ is used and if 

web length L is taken to include the curved portion of the 

web (but not the length f) in Fig. 1(b), then the results for m; = 


in reference 1 check Schuh’s results.” 
shows little effect on the 


The results! for m; = 
' check the value ¢/(1 + @), and m, 
maximum stress for the cases of g = 5.0 and 15.5. For g = 1.55, 
the variation of the stress with m; is similar to Fig. 2. For one 
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Maximuin web stresses for various joint thermal re 
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the results check Fig. 2 for 0 < » 1.0 but change to the 
@ = Scurve at m; = 

It should be noted that in many cases the variation in the 
parameters may be 


maximum stress produced by the ¢ and hy h 
no larger than the variation produced by the m; parameter 
From references 6 and 7, the range of 4; for 75ST-6 plates face to 
face is approximately: (1) 50-100 B.t.u./hour ft.? F. for asbestos 
or cement bond, (2) 225-500 for zine chromate primer or as re 
ceived or as riveted, (3) 500—-1,300 for finished surface, and (4 


700-1,300 for aluminum foil sandwich 
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INTRODUCTION 


N A SHELL or guided missile that is designed for a single use, 
I the development of appreciable creep deformations during 
the useful life need not impair proper functioning. In compar 
ison with such creep strains of the order of several per cent, 
elastic strains may often be neglected The analysis may then 
be based on a relation between the stress and the rate of creep 
strain. When the load 
of the total creep developing during the lifetime of the structure 


Dividing the lifetime into a sufhi 


intensity varies with time, an estimate 


can be obtained as follows 


ciently large number of short periods, determine the creep defor 
* The results presented in this paper were obtained in the course of re 
search sponsored by the Ballistic Research Laboratories of Aberdeen Pr 


ing Ground under Contract DA-19-020-ORD-798 
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mation that would occur during each of these periods if the 
loads were maintained at their mean values for the considered 
period. To find the desired total creep deformation, super 
impose the creep deformations obtained for the various periods 
Since the relation between stress and rate of creep strain is non 
linear, this process, as a rule, requires the solution of as many 
nonlinear boundary-value problems as there are periods to be con 
sidered. It is the purpose of the present note to show that under 
certain conditions only a single nonlinear boundary-value prob 


lem needs to be solved 
RELATION BETWEEN STRESS AND RATE OF CREEP STRAIN 


As Hoff! has pointed out, the present writer’s discussion of 
the mathematical structure of stress-strain laws in the work 
hardening range? is readily applied to the relations between stress 
and rate of creep strain 

With reference to rectangular Cartesian coordinates x; 
1, 2, 3), let o;; denote the stress tensor and e;; the strain-rate 
tensor. Since the material is treated as incompressible, 


&k = VU 


where the repeated letter subscript indicates summation in the 
customary manner. 
The stress deviation s;; is defined by 
a 9) 


$ij = oi; — (1/3) ondi; 2 


where 6;; is the Kronecker delta. A complete set of independent 


invariants of the stress deviation is given by 


The first of these vanishes by the definition of s;, \ny invariant 
of the stress deviation can therefore be written as a function of 
J,and J;. The deviation ¢;; of the square of the stress deviation 
is given by 

2/3) J6i; { 


For simple tension in the x,-direction, the relation between the 
tensile stress 0); and the corresponding creep rate €,; is generally 


given in the form 


€)) = Ca," o 


where ¢ and m are constants. For an incompressible isotropic 
material, the most general relation between the tensors of stress 
and creep strain that reduces to Eq. (5) in simple tension has the 
formT 

€i; = f(J2, Js) [¢( Je, Js) siz; + h( Je, Js) ti; 6) 


In Eq. (6) each of the functions f, g, and 4 is homogeneous in the 
components of the stress deviation, the orders of the products 
fg and fh being n — 1 and n — 2, respectively. Note that the 
creep laws suggested by Odqvist® and Bailey‘ are special cases 


of Eq. (6) 
BOUNDARY-\V ALUE PROBLEM 


Consider a body that consists of a material with the creep be 
havior defined by Eq. (6). Denote the volume occupied by the 
body by V and let its surface S consist of the parts Sy; and S7 

Starting with the stress-free, undeformed state, let the veloci 
tics v; of the points of S\- be required to vanish at all times, while 
the points of Sr are subjected to given surface tractions of the 
form 7;(x, t) = 7;*(x)¢(t), where x stands for the coordinates of 
a generic point and ¢ denotes the time. 

We will show that the stresses, velocities, and displacements 


throughout the body then have the forms 


oij(x, t) = o:;*(x)d(b) (7) 
v(x, t) = v,;*(x) V(b) 8) 
u(x,t) = u;*(x)x(b) Q) 


t For a proof of this statement see reference 2 
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where y = ¢Q’ and X - {| y dt 
0 


To prove this, we need the fundamental equations of the prob 
lem. These are (a) the equation of equilibrium, o,;,; = 0, where 
the comma subscript followed by 
b) the definition of the strain 


the subscript 7 denotes dif 


ferentiation with respect to xj; 
rate in terms of the velocity, €;; = (vi,; + 7;,:) 2; (c) the condi 
tion of incompressibility, Eq. (2); (d) the relation, Eq. (6), be 
tween strain rate and stress; (e) the boundary condition on 
stress, ojjn; = 7; on Sy, where n; is the unit vector along the 
exterior normal of S;; and (f) the boundary condition on velocity 

= Jon S; 

With the exception of the stress-strain law ot Eq 6), each of 
these fundamental equations involves, in a linear and homo 
geneous manner, either only the dynamic variables o;; and 7 
or only the kinematic variables v; and e«;, \ccordingly, the 
will drop out when Eqs. (7) and (8) are 


/ 


functions (ft) and y(t) 
substituted into these equations. The stress-strain law of Eq 
6), on the other hand, is linear in the components of the strain 
rate but homogeneous of the order » in the components of the 
This law therefore furnishes the first relation of Eq. (1 

then follows 


stress 
between ¢ and y The second relation of Eq. (19 
from Eqs. (8) and (9) and the fact that 7; = Ou; Of 

Analyzing stresses and strains in incompressible, work- 
hardening plastic solids, Ilyushin® pointed out that, under stress 
boundary conditions, the proportionality between stresses and 
loads, familiar from linear elastic solids, subsisted when the 
secant shear modulus was proportional to a power of the octa 
hedral shearing stress—i.e., a power of Jo. Aside from the sub 
stitution of ‘‘rate of creep strain”’ for ‘‘strain,’’ the result estab 
lished above is a natural extension of Ilyushin’s result 

To illustrate the effect of the exponent » in the stress-strain 
law of Eq. (5) for some typical variations of the load intensity 


with time, consider the rectangular, sinusoidal, and triangular 
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vn in Fig. 1 These pulses have the same duration 
g I 





pulses i 

( peak valu The total creep deformations pro 

1K sinusoidal or triangular pulse can be obtained by 

ipt g rrection factors As or Ar to the total creep deforma 

tion produced by the rectangular pulse. Fig. 2 shows thesc 

cort factors versus the exponent m in the stress-strain law 

f Eq For n = 1, these correction factors are the ratios 

tw iverage load intensities of the sinusoidal or triangular 

pulses | that of the rectangular pulse. The sharp decrease of 

the correction factors with increasing m is due to the fact that 

the relative contribution of loads below the peak intensity to the 
ital creep deformation is small for large values of » 
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Second Approximation to the Laminar 
Boundary-Layer Flow Over a Flat Plate? 
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Visiting Research Professor, Institute for Fluid Dynamics and Applied 


Mathematics, University of Maryland, College Park, Md 
ctober 8, 1956 


the Navier 
differential Reynolds 
was the first to consider the second approximation for the 
i.e., the Blasius flow He ob 
tained the second approximation to the stream function and dis 


I’ IS WELL KNOWN that the boundary-layer theory is the first 


approximation of the asymptotic solution of 
Stokes 
\lden 


case of a semi-infinite flat plate 


equations for large Numbers 


cussed the pressure and velocity distributions as well as the local 
The second-order correction to the 


> 


skin friction on the plate 
local skin friction turned out to be proportional to x , x being 
so that the total drag could not 
Kuo? has 


technique 


the distance from the leading edge, 
be calculated by integration of the local skin friction 
this difficulty of the Lighthill 
he considered a more difficult problem of a finite flat 


be 


by use 


overcome 
\ 


(ctually, 


instead of a semi-infinite plate. However, it should 


Kuo's second approximation is lower in degree of 


plate 
noted that 
ipproximation than Alden’s because the basic equations for his 
second approximation are nothing other than the usual boundary 

layer equations; the difference from the first approximation lies 
only in adopting the pressure distribution as calculated from the 
the \lden’s equa 


second approximation is, of course, different from 


displacement thickness of boundary layer 


tions for the 


the usual boundary-layer equations. 





nfortunately, Alden’s result contains the serious defect that 
the vorticity tends to zero as a negative power of the distance 
This 
the 





should decay exponenti 
the 


from the plate, although it 


itter requirement can easily be inferred on basis of 


Oseen equations, which are certainly valid outside the boundary 
laver 


In this paper, a second approximation to the Blasius flow 


which is free from the above-mentioned difficulty will be given 
T This investigation was carried out under Contract Nonr-595(08) of the 
Office of Naval Research, whose support is gratefully acknowledged 


On leave of absence from the University of Tokyo, Japan 
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\lso, the total skin friction acting on the of the plate of 


length x from the leading edge is calculated by us 


portion 


f the theorem 


of momentum, thus avoiding the difficulty of nonintegrabl 
singularity at the leading edge rhe procedur mploved is 


similar to that employed in the investigation of the asymptotic 
behavior of viscous fluid flow past a finite body 
\lden’s 


scheme of solution of the Navier-Stokes equations in terms of 


Recently, Bell raised a question about the validity of 


expansion In inverse powers of the lo il Reynolds Number, 
which has been suggested by von Karman Bell concluded 
that Alden’'s scheme could not lead to a solution because at the 
fourth approximation it unavoidably violates the condition at 


the outer edge of the boundary laver rhe author believes that 

this difficulty is intimately related to the above-mentioned 

slow decay of vorticity as a function of the distance from the 

plate in Alden’s solution \lso, it is believed that the original 
} 


pe success 


von Karman scheme of successive approximation can 
fully followed if appropriate consideration is given to the behavior 
of vorticity distribution in the boundary layer as shown below 


Let a viscous incompressible fluid of kinematic viscosity v be 


flowing with velocity Ul in the positive direction of the x-axis 
and a semi-infinite flat plate be placed along the positive x-axis 
with its leading edge at the origin 
In terms of the parabolic coordinates £, n de fined by 
vy = 4 
where c=x-f Y, Cc =€é- n 2 


the Navier-Stokes equations can be expressed in the form 


&* + n*) AAW — 4/5(0/0F) + n(0/O0n 1; Av = 
1/v) ;2[E(OW On n(OW, O£)} AW 
gf + 9? OW OF) (0 On OW On) (0 OF) AV, 
where A O? 0&2) + (07/0n*), and W is the stream functior 
such that the velocity components are given by 
u = OW/dy, v= —OW/da i 
It will be assumed that V can be expressed in the form 
v = ofEFin) too tut... } 
as ¢-> , Where y = ("¢ and log é is considered to be 


order terms 


Phe 


Alden tacitly assumed that even 


yo 


sent from the asymptotic expression boundary « 


are that on the plate 


and, as the edge of the boundary laver is approached, the vor 


ticity should vanish exponentially, 


— () exponentially as 7 


wo = WV Ly, 


since the boundary layer goes over smoothly to the outside irro 
tational main flow [This latter 


Alden’s and Bell’s work in favor of the 


condition has been ignored in 


Intuitive condition 


more 


of smooth joining to the uniform flow u— —ilasn— 
It should be borne in mind that the expansion in the form of 
Eq. (5) is only valid for the boundary layer where » 0 (1) but 
not for the main flow where £ and » are of the same order of 
magnitude In our case, it is also required that 

the fl ee is syvnin c tr ul u tn ré L €¢ t to tite IX Ss 


Substituting Eq. (5) into Eq. (3) and comparing terms of the 


quations for 


F is found 


same order of magnitude, we can set up a system of 
F, Yo, ¥ 


to be the Blasius function such that 


\s should be exp cted from the outset, 


F'’ + FF’ = 0 
F(0) = F(0) =0 F’'(0) = al 7 
F—+2)—s) as n— 
a = 1.32824 8 = 0.8604 


By considering the conditions of Eqs. (6), (7 


found that 
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where C,; = 0.0164, b = 1.454 
Gn «tall F | F’H/ F2F"’) dn Thus, the constant } characte rizing the friction near the leading 
J0 ~e edge is about 10 per cent higher than the Blasius value a, eo; 
7 I 0 nF’ — F)H/F?F’'|dy trary to Carrier and Lin’s conclusion that b = a 
n . 
Hin) = J (nF F)? $32} F dn + ye 
0 12) REFERENCES 
A/2) |FF 1/2) F’? 
e Aiden, H. I Second Ap ximation to the Laminar Bound La 
A - | \(nF’ F)? $32} F dn = 3.3199 Flou er a Flat Plate, J. Math. & Phys., Vol. 27, No. 2, pp. 91-104, Jy 
JV 1948 
G = —2.2975 Kuo, Y. H Fi m ) ible V Fluid pa | 
; Plate at Moderate Vumbe J. Math. & Phy Vol. 32, N 2 
and C; is an arbitrary constant 83-100, July, 1953 
Next, by use of the technique developed in the previous paper, ‘Imai, I., On the Asymptotic Behaviou Viscous Fluid Flow at a Gr 
the irrotational flow induced by the vorticity distribution in the Dista dbo ; R to Fi ore 
: a ’ Proc. Roy. Soc. London A, Vol. 208, pp. 487-516, 1951 
boundary layer can be calculated The result is that the stream ‘Bell. S.. Stud Risdinve tavay Wise! Sclnicwws wed -Aties Meth 
function can be expressed as the imaginary part of an analytic Boundary Layer Correction, University of California Institute of Eng 
function (complex velocity potential) /(¢) : neering Research Tech. Report, No. 20-106, December, 1955 
von Karman, Th., The Engtne Grappl with Nonlinear Problen Bu 
vw = Im f(¢) 13 Am. Math. Soc., Vol. 46, pp. 615-683 (especially 655-659), 1940 
Carrier, G. F., and Lin, C. C., On the Nature the Boundary La ne 
where the Leading Edge of a Flat Plate, Quart. Appl. Math., Vol. 6, No. 1, pp. f 
\ 9 68, April, 1948 
re = vié a1 P¢ C/¢ 
ipA [(log ¢/1)/¢] + O(1/¢?)} 14 
with C = 28C, + G(~) 15 + 


The local skin friction can be calculated as 


mH = pw() = (1/4) pl oy R! 2) x 
fy 1/R) [¢ A/2) log R' 7}! (16 
where R = lx/v. The first term is the well-known Blasius for 


mula. It is obvious that the ‘total skin friction’? cannot be 


calculated by the formula 


in view of the singularity at x = 0 of the term R~* 7 = O(a 
This difficulty can be avoided by considering the balance of 
momentum of the fluid contained in a large circle of radius 4 
with center at the origin. For this purpose it is only necessary 
to know the behavior of the flow on that circle, which is exactly 
to the order required) given by Eq. (14) (in the main flow) and 
by Eq. (5) (in the boundary layer). This procedure has also 


been used in the previous paper The result is 


Cr = (a/R) + (482/R) + (a/R?) (A /2) X 
l + log R’’*) — CG} + 0(1/R . 
= (1.328/R'/*) + (2.326/R) 1.328/R®’*) x - 

} 1.65901 log Ri?) 4 Cif + 0C1/R? 


where C; = D 
It is interesting to note that the lowest order correction to the 
Blasius formula is 78?/R, although an intuitive argument based 


1/2)pl?-+2x] is the total skin-friction coefficient 


on the formula for the local skin friction [Eq. (16)] would suggest 
= Also, it is remarkable that for the 


calculation of the correction term it is only necessary to know 


it to be of order O(R 


the Blasius solution if the momentum theorem is used instead of 
the integration of the local skin friction. Further, it may be 
interesting to compare the result with Kuo’s formula for th 
total skin-friction coefficient of a finite plate :* 


C, = (1.328/R'/*) + (4.12/R) (18) 


Finally, it should be mentioned that the constant C, appearing 
in Eqs. (11), (15), (16), and (17) remains undetermined in the 
above analysis. It seems to the author that this, as well as other 
constants which would appear in the higher order approximations 
Yo, vz, ..., could only be determined by taking into account the 
flow condition near the leading edge, which can be investigated 
on the basis of the Stokes-type approximations, as done by 
Carrier and Lin.* In fact, it can be shown that for small R 


Cy = (b/R'/*) — (rb2/256)R + 0(R*/?) (19 


where > is an undetermined constant. If one tentatively adopts 
the condition that C; and dC,;/dx, and hence total and local skin 
friction, should be continuous at R = 1, C; and #) can be deter 


mined as follows: 











Three-Dimensional Vortex Pattern Behind a 
Circular Cyliinder* 


Francis R. Hama 

Assistant Research Professor, Institute for Fluid Dynamics and 
Applied Mathematics, University of Maryland, College Park, Md 

October 15, 1956 


a OF THE WAKE behind circular cylinders is 
vestigated in detail by Roshko by means of a _ hot-wir 
anemometer.' The Reynolds Number range of periodic vort 


shedding is divided into the following three subranges: 


Stable range 4() lid/y 150 
Transition range 150 Ujd/v < 300 
Irregular range 300 lid/y 


is the velocity of the incoming flow, d the cylinder 
Whereas regular lami 


where Ll 
diameter, and vy the kinematic viscosity 
nar vortex streets are formed in the stable range and no turbulent 
motion develops, turbulence is initiated by laminar-turbulent 
transition in the free shear layer separated from the cylinder 
This transition begins to appear in the transition range 

On the other hand, the present writer? has investigated th 
flow phenomena associated with transition especially behind 
trip wire attached to a flat plate, finding that a two-dimensional 
discrete vortex imbedded in shear flows has a strong tendency t 
form three-dimensional vortex loops and that the formation an¢ 
development of the three-dimensional vortex loop are directly 
responsible to the creation of a turbulence spot 

In order to examine whether or not the same principle should 
apply to the transition behind a circular cylinder also, the three 
dimensional vortex pattern behind a circular cylinder is observed 
in a towing tank, 2 ft. wide, 2 ft. deep, and 12 ft. long. The circu 
lar cylinder is a 1 /4-in. outside diameter brass tube, approximatel 
20 in. long, and is submerged horizontally at the mid-depth of th 
water. Red dye is very slowly leaked through a row of small 
holes separated approximately 1/24 in. and aligned parallel to th 
cylinder axis, the diameter of the holes being 1/32 in. The row 
of holes is located approximately 120° measured from the forward 
stagnation point so that only one side of the vortex streets is ob 
served from above the plane of the vortex street 

Result of the observation is presented typically in Figs. 1 
Fig. 1 is taken at the lowest Reynolds Number among the pictures 
and stable range. The vortices are shed regularly without show 
ing much three dimensionality except for a slight long-wave de 


* This research was supported by the U.S. Air Force through the Office 
of Scientific Research, ARDC, under Contract AF 18(600)-1014 
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d/y = 117 


Vortices behind a circular cylinder ( / 


which is rather attributable to persistent large-scale 


Although lower Reynolds Numbers 


formation 
disturbance in still water 
than this could not be obtained in our towing tank, an appré 
iable long-wave three dimensionality at low Reynolds Numbers 
occurs only when the observation is made after a relatively short 
waiting period between consecutive runs evidently due to per 
otherwise reasonably good two 


sisting large-scale disturbances; 


dimensionality such as Fig. 1 is consistently obtained. We fear 
that the long-wave three dimensionality at very low Reynolds 


Numbers such as 8 reported by Roshko is due to the lack of uni 
formity of the wind-tunnel flow 
When the Reynolds Number is increased slightly larger than 


1), a three dimensionality —i. the transverse waviness 
J 
she 
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RS’ 


FORUM ey, 
begins to appear fairly randomly and quite hesitantly At tl 
Reynolds Number of around 200, the three dimensionality ap 
pears violently as shown in Fig. 2 \s the Reynolds Number 

furthermore increased, the three dimensionality appears mor 


regularly and periodically in the longitudinal as well as transvers 


directions The transverse wave length seems to decrease as the 


Reynolds Number increases In this transition range, however, 
the three dimensionality does not develop further but only rotates 
iround the original vortex axis. Consequently, the local transi 
tion from laminar to turbulent motion does not seem to occur in 
Nevertheless, a hot-wire trace should 
as observed by Roshko because 


transition occurs fat 


the shear-flow region 
indicate somewhat of irregularity 
of the 


downstream of the circular cylinder due perhaps to the interac 


three-dimensional deformation rhe 


consecutive deformed vortex lines and the two 


and to the effect in 
ments, moreover, the flow turbulence would promote these inter 


tions between 


vortex streets sid¢ wind-tunnel experi 


actions 


Flow phenomena in the irregular range, in which the Reynolds 
Number exceeds approximately 300, look much the same as those 
wire on a flat plat In Number 


shed by the circular cylinder becomes wavy 


behind a trip this Reynolds 


range, the vortex line 


almost immediately, develops further apparently in the she 


laver, and breaks down into turbulent motion— a typical picture 


being Fig. 3 


The of the three dimensionality and its important 


in the creation of the 


appearance 
role local transition into turbulent motion 
are thus confirmed in the wake 
rhese processes are therefore considered to be the guiding prin 


behind a circular cylinder also 


ciple of transition phenomena not only in boundary layers and 


et boundaries as reported elsewhere? but in wakes as well 


In connection with the above observation, it may be of interest 
to mention the destabilization effect of a thin flat plate mounted 


wake his inter 


ference device is also investigated by Roshko,! and the 


behind a evlinder in the center plane of the 
effective 


ness in stopping the periodic shedding is attributed to the inter 


ruption of the communication of the two free vortex layers 


It is true that the separation plate interrupts the communication, 


which occurs always three dimensionally, and stops a large-scale 

periodic shedding his constitutes, however, only one side of 

the effectiveness. The other and important side of the story is 
ad 
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Fic. 3. Vortices behind a circular cylinder (U,\d/vy = 318). 








POURS OF THE AERONAUTUT 


of the shear lay 
three 


ers along the flat plate, the forma 


because 


the 


that, 


tion of dimensional transverse waviness and th 


breakdown into turbulent motion are promoted and that ther 


fore, when the flow leaves the combined system of evlind: 


no p 





plate, periodic vortices are entirely destroved, a1 period 


shedding is experienced. One other feature of the effect of the 
should ntioned—that 1s, that the sep 


vnolds Num 


is extended to approximately 


separation plat ilso be m 


iration plate stabilizes the periodic vortices at low R 


\s a result, the stable rangi 


bers 


20 with a tested length of the plate, whereas the transition rang 
is extremely squeezed, and an entirely turbulent wake ts obser 
itso 
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GENERALLY ACCEPTED VIEW of boundary-laver transition 


’ i SHE 
is 


control transition if free-stream turbulence is somewhat less than 


that the Tollmien-Schlichting oscillations precede and 


0.1 per cent, while the Taylor theory, based on the concept 
of momentary separation, is only valid for free-stream turbu 
lence greater than about 0.2 per cent (c.f., Dryden rhus, all 


] 





laminar boundary 


work on 


directed toward 


of the more recent theoretical 


instability has been the improvement and ex 


tension of the Tollmien-Schlichting theory However, there is 


i part of the transition problem that seems ideally suited for an 


attack patterned after Taylor's method —-namely, transition 


caused by roughness rhis note presents results for the case of 
single roughness elements on flat plates in subsonic, low-turbu 
lence streams 


Phe 


summarized a great amount of experimental data 


present subject has been considered by Dryden,? who 


ind discussed 


several intuitively plausible parameters for correlating data 
Most of the experimental information used by Dryden was 
taken from Tani and Hama,* and some newer data have been 
published by Klebanoff, Schubauer, and Tidstrom The work 
to be described assumed model configurations and flow condi 
tions typified by these references 

In 1936, Taylor® published a theoretical discussion of the rela 


tion between scale of turbulence and the critical resistance of 


spheres. For isotropic stream turbulence produced by grids in 
subsonic wind tunnels, he was able to demonstrate that the Pohl 
suld be 


transition 


hausen parameter, VY A”, written for turbulent motion, c 


considered to attain some constant value defining 


This led him to the formulation of satisfactory parameters de 
scribing boundary-layer transition on smooth plates and spheres 
Considerable supplemental information was available to aid 
Taylor's study of transition due to isotropic turbulence produced 
by grids, but the published data on transition produced by rough 
ness do not include similar details of the flow It has been neces 


sary to cover these gaps by numerous assumptions which are 
now described 


The form of the Pohlhausen parameter to be used is 
Vr"? = 6° ul) (Op'/Ox l 
where 


6 = 


boundary-layer thickness at x 


mw = dynamic viscosity 


ICAL 
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= free-stream velocit 
= fluctuating pressure caused b scillations in 
laver 
= distance from leading edge of plat 
It is evident that ability to ¢ ite the fluctuating 
gradient due to oscil s iste ss t l 
iry layer represents the main obstacle to calculation of 
of the Pohlhauser paramete ait g insiti An 
approaci will be s \ 
Ass le that 
2) J 
= fluctuating component of velocit 
mass density 
in unknown characteristic lengt 
ilk eans “proportion lt 
Now it is stipulated that this ¢ vat S restr 
where transition is predominantly influenced by the 
element 1.¢ we egiect tree-str eyeat ¢ 
disturbing factors hus, we 1 
whe 
free-stream velocit 
location of transitiot 
= location of roughness eleme: 
k = height of roughness element 
= boundary-laver thickness corresponding t 
Phe terms in parentheses in Eq. (3 resent the lox 
size of the roughness element rhe lengt IS USE 
because Eq. (1) must be evaluated at the location of t 
Po complete Eq. (2), the lengt ust be consi 
SEC s reasonable that the lengt vetWee»l SUCCESSIVE 
shed by the roughness element should be important 
to the fluctuating pressure gradients. If this length is tak 
is s, then u’ and s should be so related that we may rewrite E 
3)as 
nu’? - Pp t vi, R/O 
Eq. (4) may be substituted into Eq. (2) s! ing that the S 
of these assumptions is 
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points plotted are for the condition 2 where Xo Corre 


: | DATA OF TAN AND HAMA sponds to k = 0, the case of no roughness. Eq. (9) does 
AERO. RES INST, TOKYO predict when x; Z x, nor does it tell when 4 Phese are 


| . ? , 
' P | SYM K~mm Xe~METERS the limiting cases not considered in the analysis, and Eq. (9 
. C Ss ¢ k~ ERS 
rm - | ° 0.50 050 probably becomes less accurate as the bounds are approached 
| oy O73 0.50 However, references 2 and 4 contain information permitting t 
1 \ o 0.73 0.80 ; ; ee 
e = te 099 050 ipproximate prediction of these limits, and they should be est 
2 \ ° 0.99 0.80 mated before using Eq. (9) for engineering 
ss | < O99 Oo 
| Fd 0.99 25 , aie 
as | 5 099 140 REFERENCES 
j & 49 0.50 ) 1H I \ 
4 49 0.80 : 
9 4 49 “ Flow,’’ A 1 ul I } 1 I 
4 x 49 40 Karman, V« pp. S-lt Academ Pre Ir New York 448 
Dryden, H R ) f 
Q S \ 20. No.7 177-48 ) 
© 12 » I I.,and Hama, F.R.,S s 


Rey X 10 
wt 
qe 
* 
- 
I 
-. 
. 


3 X Taylor, G. I., Sta t | 

THEORY neerhunielege Si Nth alee ia . 
& 1/Xs I= - 9 

— Rey (oN st)” const 7-3 

jn a (Xe Xa7 = 780 


An Approximate Method Applied to the 
Solution of the Problem of Vibrating 
Rectangular Plates 














tr sit } 0 7 " Milom y otanisic 
2 5 0 = @) 25 e a5 
O 5 0 5 fe 25 30 35 , 
K ix / A oc ate P fe .o c nw ee 
/ \ = 1/ 
LS es x 104 Purdue Unive Lata nd 
Xk \ Xk ; 
in reg October 1, 195 
is tak Fic. 2 
r ' 
vrite | i SUMMARY 
Op’, Ox 2 vk, R/ 7) : 
Ordinary product tion and Galerkin method are i 
Substituting the above result in Eq 1) vields the expression the natural frequency of a rectangular s tl 
boundary conditions a) all edges fixed woa xe 
the resu j V XX"? b- ul pl 2/k) folxs/xx, 2/6 6 and the other two edges free The plates n i ime 
»bey the classical small deflection theory, and the influence of rotatory 
is to be evaluated at 2 ind 6; at xXx, so Eq 6) becomes tia and shear forces are neglected in the solution R ts are ta ited 
Vv k Vi Xk, R vA ‘ SYMBOLS 
4 ° 
mplete the right side of Eq. (7) and find if the result is rectangular coordinate 
rm constant at transition, it is necessary to determine the lateral dimensions at I 
rs 
1 , ; trans sal deflection tl 
us inside brackets. Inspection of published data, oes i 
. : thickne f ate 1 
rst for vi = const., indicates that , } 1)] — 
aie ‘ 
R R vr “| = (x:/R Rk VX - I n's ra 
mass density of ite materia »s € 
This1 be rewritten as angular frequency of conservative syste " 
angular frequency of the nonconservative system, rad 
Ket rR 2 ~ const ‘ VXy 2 const \ ph D) « value defined by | 
A \ eigenva t 
rat = x 
i = kinematic viscosity as _ - — 
Having obtained this information, it is easy to establish fo(x;/a ; Ma ———- > X 
ete the parameter We find that 
. > / 
V A R R/x uy /X - Ss / 
- ae 
Re k/x Ve /X “= = const Q j 
j 
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single roughness elements on flat plates in low-speed streams 





I validity of tl derivation cannot be checked without / 
experimental measurements. However, the result ex- / 
ssed by Eq. (9) is readily verified by available data is demon- / 
strated by Figs. 1 and 2. From these figures, it is found that / 


¢ em ° ony ° . » » r 
—o mstant of Eq. (9) is approximately 780 for both sets of data ; 
Here it should be noted that the data points have been read a _ 


) 80 from small figures in the referenced publications, so it is possible y 





ill and probably consistent errors may be present 
I ; : ; : ; 
' rhese might make a slight change in the value of the constant, 
ugh this is unlikely in view of the experimental scatter Phe Fic. 1. Coordinate system of the plate 





JOURNAL OF THE 


positive integers 


coefficient used in series representation of deflection 
a function of x and y alone 


a function of x alone 


a function of y alone 
time-depending function 
characteristic function of the vibrating beam 
damping coefficient 
a coefficient defined by Eq. (10 
é, Kronecker delta 


a 
Derivatives with respect to x and y are indicated by subscript x, xa 


and, with respect to time, by ¢, ¢/ 


SOLUTION OF THE PROBLEM 


——* ro the classical two-dimensional theory of flexural 
motion of a vibrating plate, assuming that the damping 
forces are proportional to the velocity, neglecting rotatory in 
ertia and shear forces, the motion of the plate is governed by the 


following partial differential equation: 
Atw(x, y, t) + (ph/D)wilsx, y, t) = —(k/D)wix, y, t) (1) 


Using an ordinary product solution, such as 


Z 
w(x, y,t) = 2. onl x, y(t) 
m= 0 
the homogeneous equation 


Viw(x, y, 1) + (ph/D)wilx 


becomes V 4on(x, y) Agn( xX, ¥) 
Then Vig» (x, ¥) = (ph/D)w no, 
Hence, Eqs. (1), (2), and (5) lead to 
W(t) + (R/ph) Vi, lt 
the solution of which is 
W,,(t) = A,e~@ cos (Q,,¢ 
Eqs. (6) and (7) show that 
2a (k/ph) = 0 (S) 
a® + (k/ph)a m* (9) 
Therefore, = k/2ph (10) 
[wm? — (k/2ph)2]'/* (11) 
Hence the solution is determined as soon as the angular frequency 


of vibration of the conservative system w,, is found. In other 
words, the eigenvalue of Eq. (4) must be determined 


A most convenient method is that of Galerkin. Using 


] j 
d(x, y) = rh > ye eS 
) Ea 1 


where Y and Y are of the form 
n(&) = cosh p,(&//) cos p,(é/1) — 


y, {sinh p,(é//) sin p,(£//) 


the Galerkin form of the stated problem becomes 


i j a fb 
} a 3% taf | [V4(X,¥,) — MX-V)|Xp¥ydxdy =0 (14) 
1 s=1 0 0 


RERONAUTICAL SCIENCES 


TABLE 1] 
lhe Eigenfrequencies, Al/?, of the Plate for Different Shape 
Ratios and k = 0 


Type of Plate Support 


All edges fixed 


Pwo adjacent edges fixed and 
the two other edges free 


Two cases will be considered in this paper—namely, (a) a 
plate with all edges fixed, and (b) a plate with two adjacent edges 
fixed and the other two edges free. The characteristic functions, 
Eq. (13), are orthonormal and satisfy the boundary conditions— 


namely, 
(a) n (0) 
(b) 70 


Then Eqs. (13) and (15) lead to 


cosh By COS pb 

= (cosh M COS pb (sinh MK sin lu 
while Eqs. (13) and (16) lead to 

cosh cM, COS p 

Yr = (cosh p, + cos z,)/(sinh pw, + sin pg 


With the abbreviated notations 


and, since the characteristic functions are orthonormal, it follows 
that 


[ fan 


where 
l for p = randgq = 


0 for either p # rorg # s 
ES = (b/a ur? + (a/b)3y.4 + 2H,,-H.,, for p=r ind g=s 
E,, = 2(a/b)H,,H,, for either p ¥ rorq # s 


Finally, since = phw,,2/D 


or wm? = AD/pha*b 


then 2, = [(AD/phasb) (k/2ph)2|'/* 


The results for uw, and y, are given in any standard reference of 
vibration. 

The coefficients Ez can easily be determined for various 
shape ratios and for each boundary condition. Using only two 
terms in the series representation, the eigenfrequencies, Al ~| for 
each case under consideration and for k = 0 are shown in Table 1. 
The results show that the mode number is irregular with in- 


creasing shape ratio 





ent edges 
functions, 


iditions— 


it follows 


19) 


erence of 


r various 
only two 
Al’, for 
Table 1. 


with in- 





